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1. 


INTRODUCTION 


This  report  documents  the  work  performed  under  a  Phase  II  SBIR  study  (of  2  year 
duration)  sponsored  by  AFOSR.  The  objective  of  the  overall  project  (Phases  I  and  II) 
is  to  develop  an  experimental/computational  approach  to  model  precombustion 
chemistry  and  transport  behavior  of  hydrocarbons  under  supercritical  conditions. 

The  model  development  and  computations  were  performed  at  CFDRC  and  the 
parallel  experimental  effort  was  conducted  by  Prof.  L.D.  Chen  at  the  University  of 
Iowa  (under  a  subcontract  from  CFDRC).  The  purpose  of  this  report  is  to  provide  a 
recapitulation  of  the  project  objectives,  overview  of  the  technical  approach  and  a 
detailed  description  of  the  progress  made  during  the  study. 

1.1  Significance  of  the  Topic 

Modern  military  aircraft  employ  fuel  as  the  primary  heat  sink  medium  for  heat 
loads  arising  from  sources  such  as  the  engine,  the  avionics,  the  environmental 
control  system,  and  the  air  frame.  It  is  expected  that  advances  in  engine  technology 
for  high  performance  aircraft  are  going  to  further  increase  the  heat  loads  on  the  fuel. 
The  increase  in  fuel  temperature  due  to  heating  has  several  consequences.  Some  of 
these  are: 

a.  Current  hydrocarbon  fuels  have  a  practical  temperature  limit  beyond  which 
the  fuel  will  chemically  decompose  to  form  gums  and  insolubles  that  stick  to 
fuel  lines  and  valves.  The  fouling  of  heat  exchangers  and  fuel  nozzles  can 
cause  serious  problems  in  the  fuel  system  [1,2]. 

b.  Pressures  in  current  fuel  systems  are  generally  above  the  critical  pressure  of 
the  fuel.  Large  heat  loads  can  cause  the  fuel  temperature  to  increase  beyond 
the  critical  temperature  of  the  fuel.  This  necessitates  the  operation  of  the  fuel 
in  the  supercritical  regime.  The  supercritical  properties  of  substances  are  not 
well  understood.  Thermophysical  properties  such  as  density,  conductivity, 
viscosity  and  specific  heat  undergo  large  variations  near  the  critical  point. 
Furthermore,  the  understanding  of  fuel  stability  under  supercritical 
conditions  is  almost  nonexistent  [3]. 
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Therefore,  efficient  heat  management  in  high  performance  aircraft  will  require  good 
understanding  of  the  transport  and  chemical  properties  of  the  fuel  under 
supercritical  conditions. 

1.2  The  Thermodynamic  Critical  State 

The  state  of  a  simple,  compressible  substance  can  be  determined  by  specifying  two 
state  variables-pressure  and  temperature,  for  example.  For  any  such  substance,  the 
choice  of  these  two  variables  determines  its  state  (e.g.,  solid,  liquid,  or  gas).  The  solid 
state  is  distinct  from  the  two  fluid  states  in  that  the  substance  maintains  a  constant 
shape  (does  not  deform  continuously  when  subjected  to  a  shear  stress).  The 
transition  from  a  solid  to  a  fluid  is  marked  by  the  loss  of  this  property.  At  low 
pressures,  the  transition  from  one  fluid  state  to  the  other  is  marked  by  a 
discontinuous  change  in  properties  at  the  boiling  point.  At  high  pressures, 
however,  there  is  no  distinct  transition  from  a  liquid  to  a  gas,  and  properties  vary 
continuously  from  one  to  the  other.  The  relationship  between  pressure, 
temperature,  and  state  for  a  typical  substance  (that  contracts  on  freezing)  is  shown 
graphically  in  Figure  1-1. 


Figure  1-1.  Phase  Diagram 
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At  low  pressures,  increasing  the  temperature  of  a  liquid  will  cause  it  to  boil  at  a 
certain  temperature,  known  as  the  saturation  temperature  (Tsat).  This  process 
involves  a  discontinuous  change  of  state  from  a  liquid  to  a  gas,  with  the 
accompanying  discontinuous  change  in  properties.  As  the  pressure  is  increased,  the 
boiling  point  increases  also,  and  the  difference  between  the  properties  of  the  liquid 
and  gas  becomes  smaller.  At  a  certain  pressure,  the  boiling  point  reaches  a 
maximum,  and  the  difference  between  the  properties  of  the  liquid  and  gas  becomes 
zero.  At  this  pressure  and  temperature,  the  substance  is  in  the  critical  state.  At 
pressures  higher  than  the  critical  pressure,  boiling  does  not  exist,  and  increasing 
temperature  produces  rapid,  but  continuous,  property  changes. 

When  dealing  with  near-critical  fluids,  state  variables  are  commonly  stated  in 
relation  to  the  critical  point.  The  ratio  of  a  property  at  a  given  state  to  the  same 
property  at  the  critical  state  is  known  as  the  reduced  property.  Also,  supercritical 
means  having  a  pressure  greater  that  the  critical  pressure  and  a  temperature  greater 
than  the  critical  temperature,  or  a  reduced  pressure  and  temperature  greater  than 
one,  and  subcritical  means  having  a  pressure  less  than  the  critical  pressure  and  a 
temperature  less  than  the  critical  temperature,  or  a  reduced  pressure  and 
temperature  less  than  one.  Near-critical  is  taken  to  mean  having  a  reduced  pressure 
and  temperature  slightly  greater  or  slightly  less  than  one. 

1.3  Property  Variations  near  the  Critical  State 

The  most  significant  aspect  of  dealing  with  near-critical  fluids  is  the  fact  that  the 
thermodynamic  and  transport  properties  vary  rapidly  with  temperature  [4]. 
(Properties  vary  with  pressure  also,  but  in  practical  situations,  variation  due  to 
pressure  gradients  is  typically  very  small  compared  to  variation  due  to  temperature 
gradients.)  Some  properties  can  even  become  infinite  exactly  at  the  critical  point. 
For  example,  density,  thermal  conductivity,  and  viscosity  all  decrease  with 
increasing  temperature,  decreasing  most  rapidly  near  the  critical  point. 
Experimental  work  has  revealed  an  anomaly  in  the  thermal  conductivity,  however. 
Investigations  involving  many  different  fluids  have  revealed  a  large  spike  in  the 
thermal  conductivity  near  the  critical  point  [5].  The  specific  heat  at  constant 
pressure  also  shows  a  large  spike  at  near-critical  conditions,  becoming  infinite 
exactly  at  the  critical  point.  This  infinity  can  be  compared  to  the  infinite  specific  heat 
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(slope  of  enthalpy  vs.  temperature)  at  subcritical  boiling  points.  The  supercritical 
temperature  at  which  the  specific  heat  is  a  maximum  is  known  as  the  pseudocritical 
temperature  (Tpc),  similar  to  the  subcritical  boiling  point.  In  addition  to  large 
property  variations,  another  important  aspect  of  near-critical  fluids  is  that  surface 
tension  is  very  small  at  high  subcritical  pressures,  and  absent  at  supercritical 
pressures. 

A  central  property  of  concern  in  fluid  flow  problems  is  specific  enthalpy,  defined  as 
specific  internal  energy  plus  the  product  of  pressure  and  specific  volume.  This  is 
important  because  it  accounts  for  flow  work.  If  all  other  variations  in  energy  (e.g. 
kinetic,  gravitational,  potential)  are  negligible,  the  variation  in  enthalpy  of  a  fluid 
will  simply  equal  the  heat  transfer  to  the  fluid.  In  this  case,  the  enthalpy  at  any 
point  in  an  internal  flow  process  will  equal  the  enthalpy  at  the  inlet  plus  the  heat 
added  between  the  inlet  and  the  point  in  question. 

1.4  Studies  of  Heat  Transfer  near  the  Critical  Point 

The  majority  of  investigations  of  near-critical  heat  transfer  phenomena  have  been 
concerned  with  flow  through  heated  pipes.  This  type  of  flow  is  not  necessarily 
forced  convection,  since  it  has  been  determined  that  buoyancy  effects  can  play  a  large 
role  in  near-critical  heat  transfer  (see  e.g.  [6]).  Thus  data  on  near-critical  heat  transfer 
are  more  properly  separated  into  two  categories:  purely  forced  convection  and  mixed 
convection.  For  large  flow  rates  in  small  pipes,  buoyancy  forces  are  negligible,  and 
forced  convection  occurs.  For  larger  pipes  or  lower  flow  rates,  buoyancy  forces  begin 
to  have  a  significant  effect,  and  mixed  convection  occurs. 

Hall  and  Jackson  [6]  cite  several  proposed  criteria  for  determining  whether  or  not 
buoyancy  effects  are  important.  They  also  note  that  all  of  the  criteria  possess  the 
form  of  a  ratio  of  Grashof  number  over  Reynolds  number.  The  Grashof  number 
represents  a  comparison  of  buoyancy  forces  to  viscous  forces,  and  the  Reynolds 
number  represents  a  comparison  of  inertia  forces  to  viscous  forces.  Thus  the  ratio 
represents  a  comparison  of  buoyancy  forces  to  inertia  forces.  In  a  more  recent 
publication,  Polyakov  [7]  presents  a  more  elaborate  criterion,  which  compares  the 
Grashof  number,  Reynolds  number,  and  a  number  which  accounts  for  thermal 
acceleration  effects. 
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When  buoyancy  forces  are  negligible,  the  flow  does  not  depend  on  the  orientation 
with  respect  to  gravity.  In  this  relatively  simple  case,  the  heat  transfer  characteristics 
are  perhaps  easiest  to  understand.  Heat  transfer  at  near-critical  conditions  has  been 
regarded  as  a  region  where  boiling  and  convection  merge  [4].  Some  aspects  of  near- 
critical  forced  convection  may  be  understood  by  examining  what  happens  as  the 
critical  pressure  is  approached  from  either  side. 

At  highly  subcritical  pressure,  the  flowing  fluid  will  change  from  liquid  to  vapor  by 
passing  through  the  well-known  regimes  of  nucleate  boiling,  annular  flow,  and  will 
eventually  reach  dryout.  This  is  illustrated  in  Figure  1-2,  which  shows  a  simplified 
sketch  of  fluid  and  tube  wall  temperatures  for  a  typical  flow  boiling  process  (constant 
heat  flux,  vertical  flow)  well  below  the  critical  point  [8].  In  the  initial  single-phase 
region,  the  two  temperatures  rise  steadily  with  each  other.  As  the  liquid  begins  to 
boil,  however,  the  wall  temperature  begins  to  decrease,  while  the  fluid  temperature 
remains  constant  at  the  boiling  point.  The  wall  temperature  continues  to  decrease 
through  the  nucleate  boiling  and  annular  flow  regimes.  When  annular  flow  dryout 
occurs,  the  wall  temperature  jumps  sharply.  When  all  of  the  liquid  has  evaporated, 
and  the  flow  again  enters  single  phase  convection,  the  two  temperatures  rise 
steadily  again. 


Figure  1-2.  Characteristics  of  Temperature  Profiles  for  a  Typical  Boiling  Process 
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Single  phase  forced  convection  data  are  usually  correlated  in  terms  of  a  heat  transfer 
coefficient.  The  validity  of  a  heat  transfer  coefficient  for  fluids  with  widely  varying 
properties  has  been  questioned,  however  [6,9].  Increasing  heat  flux  would  increase 
the  temperature  gradient,  and  thus  the  property  variation,  within  the  boundary 
layer,  affecting  the  heat  transfer.  Thus  the  heat  transfer  coefficient  becomes  a 
function  of  heat  flux.  The  concept  has  nonetheless  been  widely  used. 

One  feature  of  supercritical  convection  is  that,  in  contrast  to  subcritical  boiling,  there 
is  no  heat  transfer  crisis.  Another  is  that  the  heat  transfer  coefficient  increases  when 
the  wall  temperature  is  above  the  pseudocritical  temperature,  and  the  bulk  fluid 
temperature  is  slightly  below.  Hall  and  Jackson  [6]  demonstrate  that  for  fixed  flow 
conditions  and  low  heat  fluxes  (small  property  variations),  the  heat  transfer 
coefficient  will  depend  primarily  on  the  specific  heat  of  the  fluid.  If  the  fluid  is  near 
the  pseudocritical  temperature,  where  there  is  a  peak  in  the  specific  heat,  the  heat 
transfer  coefficient  should  have  a  peak  also.  Hall  and  Jackson  also  cite  data  showing 
that  this  is  the  case.  However,  the  peak  in  the  heat  transfer  coefficient  decreases 
with  increasing  heat  flux,  which  is  consistent  with  a  contracting  layer  of  high 
specific  heat  fluid.  Supercritical  convection  becomes  much  more  complicated  at 
higher  heat  fluxes,  when  variations  in  viscosity  and  thermal  conductivity  affect  the 
flow  and  heat  transfer  processes. 

Thus  near-critical  heat  transfer  in  the  absence  of  buoyancy  effects  seems  to  be  a  cross 
between  film  boiling  and  boundary  layer  flow  with  large  property  variations  in  the 
boundary  layer.  Heat  transfer  becomes  more  complicated  when  buoyancy  forces 
begin  to  have  a  significant  effect,  and  mixed  convection  occurs.  In  this  case,  one 
might  expect  heat  transfer  to  be  improved  in  upward  flow,  and  degraded  in 
downward  flow,  but  in  fact  the  opposite  is  true.  A  characteristic  feature  of  buoyancy 
affected  upward  flow  is  increased  wall  temperatures,  and  a  phenomenon  known  as 
the  inlet  wall  temperature  peak,  in  which  the  tube  surface  temperature  near  the 
inlet  exhibits  a  sharp  peak.  Data  also  indicate  a  continuous  improvement  in  heat 
transfer  in  buoyancy  affected  downward  flow  as  buoyancy  forces  become 
progressively  stronger  [6]. 

A  recent  study  by  Kurganov  and  Kaptilnyi  [10]  has  shed  some  light  on  these 
phenomena.  Detailed  measurements  were  made  of  velocity  and  temperature 
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profiles  across  a  vertical  tube  under  conditions  of  forced  and  mixed  convection.  For 
forced  convection  and  downward  mixed  convection,  the  flow  has  the  usual  convex 
shape.  For  upward  mixed  convection,  however,  the  flow  develops  an  M-shaped 
velocity  profile.  There  is  a  thin  layer  of  low  density,  low  viscosity  fluid  against  the 
wall,  which  moves  much  more  rapidly  than  the  core  of  higher  density  fluid.  As  the 
flow  progresses,  the  core  temperature  increases,  and  the  M-shape  smears  out.  It  is 
the  initial  development  of  this  M-shaped  profile  near  the  inlet  which  causes  the 
temperature  peak.  It  was  concluded  that  the  reduction  or  enhancement  of  turbulent 
transport  plays  a  large  role  in  buoyancy  affected  flows. 

The  effect  of  buoyancy  in  horizontal  flow  is  to  make  the  top  of  the  tube  hotter  than 
the  bottom.  Experiments  under  these  circumstances  thus  require  temperature 
measurements  around  the  perimeter  of  the  tube  to  be  useful. 

1*5  Pressure  Oscillations  Associated  with  Boiling  and  Supercritical  Heating 

Experimental  work  on  heat  transfer  to  supercritical  fluids  has  been  proceeding  for 
several  decades,  and  a  somewhat  frequent  occurrence  is  the  phenomenon  of 
thermo-acoustic  oscillations  [11].  These  oscillations  are  characterized  by  large 
amplitude  oscillations  in  fluid  pressure,  which  are  strong  enough  to  be  audible,  and 
to  sometimes  fracture  or  even  rupture  the  tubing.  Most  observed  instances  occur  in 
forced  flow,  but  they  have  been  observed  in  free  convection  also.  Pressure 
oscillations  have  been  reported  for  internal  flow  conditions  for  a  wide  range  of 
fluids,  including  hydrogen,  helium,  hydrocarbon  fuels,  and  water. 


The  conditions  under  which  these  oscillations  occur  seem  to  be  extremely  high  heat 
fluxes  and  that  the  pseudocritical  temperature  lie  between  the  bulk  fluid 
temperature  and  the  wall  temperature  of  the  pipe.  The  oscillations  have  also  been 
observed  in  a  few  instances  of  subcritical  boiling  [12],  in  which  case  the  boiling  point 
fell  between  the  bulk  fluid  temperature  and  the  tube  wall  temperature.  These 
considerations  seem  to  imply  that  the  oscillations  are  associated  with  rapid  property 
(or  phase)  change  of  the  fluid  near  the  wall. 

The  observed  frequencies  vary  widely  depending  on  experimental  conditions,  from 
fractions  of  Hertz  to  thousands  of  Hertz.  Also,  most  researchers  report  improved 
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heat  transfer  accompanying  the  oscillations.  Goldmann  [9],  for  example,  reports  that 
oscillations  in  supercritical  water  were  accompanied  by  a  decrease  in  tube  wall 
temperatures,  indicating  improved  heat  transfer.  Hines  and  Wolf  [13]  also  report 
increased  heat  transfer  accompanying  pressure  oscillations. 

Several  attempts  have  been  made  to  explain  this  phenomenon.  In  a  1961  paper  [9] 
Goldmann  pictured  packets  of  dense,  liquid-like  fluid  exploding  into  vapor-like 
fluid  upon  coming  into  contact  with  the  wall,  agitating  the  liquid-like  core,  and 
causing  more  packets  to  contact  the  wall,  Malkina  et.  al.  [14]  described  a  process  in 
which  vapor  or  vapor-like  bubbles  depart  from  the  wall  and  collapse  again  in  the 
cooler  space  of  the  core.  It  seems  likely  at  this  point  that  the  oscillations  are  due  to 
resonance  of  acoustic  waves  in  the  pipes.  Some  experimental  work  such  as  that  of 
Hines  and  Wolf  [13],  indicate  that  this  is  not  the  case,  but  there  seems  to  be  a  large 
and  growing  body  of  work  indicating  that  it  is. 

Stewart  et  al.  [11],  have  written  a  paper  in  which  they  conduct  a  brief  literature 
review  concerning  pressure  oscillations  before  going  on  the  describe  their  own 
experiments  aimed  specifically  at  studying  this  phenomenon.  They  first  cite  a  study 
involving  free  convection  in  a  vertical  heated  pipe  which  was  part  of  a  closed 
circulation  loop.  By  varying  the  bulk  fluid  temperature  (and  thus  the  speed  of 
sound),  or  the  total  length  of  the  loop,  and  observing  the  resulting  frequency  of 
oscillation,  the  researchers  found  that  the  oscillations  were  consistent  with  standing 
pressure  waves  resonating  in  the  loop.  From  their  own  experiments  with 
supercritical  water,  Stewart  et  al.  conclude  that  the  pressure  oscillations  were  due  to 
longitudinal  acoustic  waves  standing  in  the  test  section.  Their  experimental  setup 
had  voids  in  the  flanges  at  both  ends  of  the  test  section,  and  the  acoustic  waves  were 
consistent  with  open-open  pipe  resonance.  The  researchers  also  calculate  that  the 
heat  transfer  coefficient  is  a  minimum  at  locations  where  the  pressure  fluctuation  is 
maximum,  and  the  coefficient  is  a  maximum  at  locations  where  the  pressure 
fluctuation  is  a  minimum.  This  may  explain  the  discrepancy  among  previous 
reports  of  improved  or  degraded  heat  transfer.  They  also  attribute  the  driving  force 
behind  the  oscillations  to  alternating  compression  and  expansion  of  the  thin  vapor¬ 
like  layer  adjacent  to  the  wall,  producing  a  net  work  output.  The  results  of  this 
work  [11]  are  very  encouraging,  but  still  do  not  explain  the  full  range  of  observed 
cases  of  thermo-acoustic  oscillations. 
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A  recent,  very  brief  report  by  Kalbaliev  and  Verdiev  [15]  describes  the  results  of  tests 
performed  with  supercritical  toluene.  Correlations  are  given  for  heat  transfer,  wall 
temperature,  and  conditions  for  the  onset  of  pressure  oscillations,  but  no  discussion 
is  given  on  the  applicability  of  these  correlations  to  other  fluids  or  flow  conditions. 

1.6  Predictive  Capabilities  for  Supercritical  Fluids/Flows 

The  efficiency  of  using  fuel  as  an  onboard  coolant  in  military  aircraft  will  depend  on 
the  heat  transfer  and  chemical  characteristics  of  the  fuel.  These  phenomena,  in 
turn,  depend  on  the  thermophysical  properties  of  the  fluid.  It  has  already  been 
emphasized  that  thermophysical  properties  undergo  steep  variations  near  the 
critical  point.  These  variations  are  expected  to  significantly  alter  local  heat  transfer 
rates.  Therefore,  conventional  heat  transfer  correlations,  based  on  the  assumption 
of  uniform  transport  properties,  will  not  be  valid  for  supercritical  environments. 
The  use  of  ideal  gas  models  for  computing  density  and  specific  heat  is  not 
appropriate  for  the  high  pressure  conditions  in  supercritical  flows.  Also,  it  is  not 
clear  how  the  mutual  interactions  between  these  properties  will  affect  the  overall 
transport  process. 

Experimentation  with  supercritical  fluids  has  provided  understanding  on  certain 
global  aspects  of  the  transport  process.  However,  event  the  most  detailed  and 
sophisticated  experiments  are  unable  to  provide  knowledge  of  local  effects  and 
interactions  in  supercritical  flows.  Also,  most  experiments  are  done  in  simple 
configurations  while  the  applications  of  these  flows  involve  complex  geometrical 
effects  such  as  those  encountered  in  heat  exchanger  coils,  valves,  etc.  However,  it  is 
possible  to  account  for  these  complex  effects  through  the  development  of  advanced 
models.  The  great  advantage  of  modeling  is  that  detailed  parametric  studies  can  be 
performed  at  a  fraction  of  the  time  and  cost  associated  with  hardware 
experimentation.  The  model  should  possess  the  following  features  in  order  to  be  of 
practical  utility  in  designing  cooling  systems  for  aircraft: 

a.  The  model  should  calculate  all  thermophysical  properties  locally  as  a 
function  of  pressure,  temperature  and  composition.  These  properties  should 
be  input  into  a  general  purpose  transport  equation  solver  in  order  to  account 
for  mutual  interactions  between  these  effects. 
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b.  The  transport  equation  solver  should  be  able  to  account  for  complex 
geometrical  effects.  Flows  through  heat  exchanger  coils,  fuel  lines,  nozzles 
and  valves  should  be  modeled  efficiently  and  accurately. 

c.  Additional  complexities  such  as  flow  turbulence  and  chemical  decomposition 
of  the  fuel  should  be  modeled. 

d.  The  supercritical  models  and  the  transport  equation  solver  should  be 
validated  against  benchmark  experiments. 

A  model  with  the  above  mentioned  features  will  be  a  valuable  tool  in  improving 
existing  designs  and  developing  new  optimal  designs  for  aircraft  cooling  systems. 

1.7  Proposed  Phase  II  Program 

The  main  objective  of  the  Phase  II  program  is  to  produce  an  experimentally 
validated  computational  tool  for  designing  heat  management  systems  for  high 
performance  aircraft.  The  following  developments  and  improvements  were 
proposed  for  the  Phase  II  study. 

a.  Computational  Model  Development:  The  supercritical  transport  models 
developed  during  Phase  I  will  be  extended  to  include  a  range  of  hydrocarbon 
fuels.  The  models  will  then  be  incorporated  into  a  multi  dimensional,  multi¬ 
block,  time  accurate  CFD  code. 

b.  LES  Model  Incorporation:  An  LES  model  will  be  used  to  address  issues  of 
turbulence  and  unsteadiness  in  supercritical  flows. 

c.  Chemistry  Model  Incorporation:  The  chemistry  models  developed  at  Wright 
Laboratory  (for  studying  thermal  stability  of  fuels)  will  be  incorporated  into 
the  code.  These  models  will  be  further  calibrated  and  refined  for  high 
pressure  conditions  based  on  data  obtained  during  this  program. 

d.  Experimental  Testing:  Benchmark  heat  transfer  experiments  will  be 
performed  at  the  University  of  Iowa  and  Wright  Laboratory.  The 
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experiments  at  the  University  of  Iowa  will  involve  the  use  of  a  simulant 
fluid  such  as  SF6  while  those  at  Wright  Lab  will  involve  experimentation 
with  hydrocarbon  fuels. 

e.  Validation  and  Parametric  Studies:  After  the  code  has  been  satisfactorily 
validated,  it  will  be  used  to  conduct  a  series  of  parametric  studies  to  assess 
effects  of  various  operating  and  boundary  conditions  on  the  transport  process. 

f.  Customization/Documentation/Packaging  of  the  Code:  The  software 
capability  developed  in  Phase  II  will  be  customized  and  packaged  for  design 
analysis  studies.  Detailed  documentation  will  be  provided  for  usage  of  the 
code  along  with  technical  descriptions  of  the  models. 

1.8  Summary  of  Accomplishments 

The  following  tasks  have  been  completed  during  the  Phase  II  study  and  are 

discussed  in  greater  detail  in  the  following  sections  of  this  report: 

1.  Detailed  validation  studies  were  performed  to  assess  the  accuracy  of  the 
transport  models  developed  in  Phase  I.  Properties  such  as  viscosity, 
conductivity,  density  and  enthalpy  were  computed  (as  a  function  of 
temperature  and  pressure)  for  propane,  carbon  dioxide,  JP-5  and  SF6  under 
sub-critical,  critical  and  supercritical  conditions.  Good  agreement  was 
obtained  with  data  over  a  range  of  pressure  and  temperature  in  the  near 
critical  regime.  The  models  for  thermophysical  properties  were  incorporated 
into  CFDRC's  general  purpose  Computational  Fluid  Dynamics  code,  CFD- 
ACE. 

2.  The  transport  model  was  applied  to  compute  flow  and  heat  transfer  in  generic 
engineering  configurations  for  which  reliable  data  was  available.  The 
following  cases  were  considered: 

(i)  Laminar  free  convection  flow  of  supercritical  carbon  dioxide  over  a  flat 
plate. 
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(ii)  Laminar  forced  convection  flow  of  supercritical  carbon  dioxide  for  a  jet 
impinging  on  a  flat  plate. 

The  predicted  heat  transfer  rates  agreed  well  with  experimental  data  and 
theoretical  results  for  the  above  cases. 

3.  The  two  layer  k-e  model  for  turbulence  was  incorporated  into  the  code.  It  was 
tested  on  a  number  of  benchmark  problems  (for  incompressible  flows)  and 
the  results  were  compared  with  those  predicted  by  the  standard  k-e  model. 
The  mean  velocity  profiles,  turbulent  kinetic  energy,  pressure  coefficients  and 
skin  friction  coefficients  predicted  by  both  models  were  in  close  agreement  for 
all  the  cases  considered. 

4.  The  low  Reynolds  number  k-e  model  was  used  successfully  to  model  heat 
transfer  and  chemistry  in  supercritical  fluids.  However,  a  few  advanced 
features  had  to  be  incorporated  into  this  model  to  account  for  steep  property 
variations  in  the  fluid.  An  advanced  correlation  was  used  to  calculate  the 
local  turbulent  Prandtl  number  to  account  for  the  "bursting"  phenomenon 
observed  in  the  sublayer.  With  this  addition,  the  low  Reynolds  number 
model  was  able  to  resolve  experimental  features  that  could  not  be  modeled 

using  the  standard  k-e  model. 

5.  The  two-layer  k-e  model  was  tested  for  heat  transfer  applications.  Its 

predictions  were  found  to  be  less  accurate  than  those  of  the  standard  k-e 
model  and  the  low  Reynolds  number  model.  This  was  mainly  due  to  the  lack 
of  an  appropriate  heat  transfer  analog  (of  the  momentum  transfer  model)  for 
the  two-layer  model. 

6.  The  models  were  used  to  simulate  laminar  and  turbulent  heat  transfer  in 
channel  flows  and  tube  flows  in  the  supercritical  regime.  The  simulations 
were  done  for  carbon  dioxide  since  a  large  amount  of  experimental  data  exists 
for  the  flow  and  heat  transfer  characteristics  of  supercritical  C02.  The 
importance  of  buoyancy  in  influencing  local  heat  transfer  rates  was 
investigated.  The  results  predicted  by  the  model  were  found  to  be  in  good 
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quantitative  agreement  with  existing  data.  A  detailed  parametric  study  was 
conducted  over  a  range  of  Reynolds  and  Grashof  numbers  to  determine  the 
effects  of  natural  convection. 

7.  Experiments  were  performed  at  the  University  of  Iowa  (under  the 
supervision  of  Prof.  L.D.  Chen)  using  SF6  as  a  surrogate  fuel.  During  the 
summer  months  (June  to  August,  1994),  experiments  were  also  performed  at 
Wright  Laboratory  to  obtain  deposition  data  on  fuels.  This  data  was  used  for 
validation  studies  of  the  computational  models  developed  during  the  Phase 
II  program. 

8.  The  Krazinski  model  for  the  auto-oxidation  of  fuel  was  incorporated  into  the 
code.  The  results  predicted  by  the  model  were  in  agreement  with  data 
obtained  by  Marteney  and  Spadaccini  as  well  as  predictions  by  Krazinski  et  al. 

9.  Advanced  multi-step  models  for  thermal  stability  of  hydrocarbon  fuels 
(developed  at  Wright  Laboratory)  were  incorporated  into  the  code.  In 
particular,  a  9-step  model  (consisting  of  6  gas  phase  and  3  surface  reactions) 
developed  at  Wright  Laboratory  was  used  in  the  simulations.  The  model 
predictions  of  wall  temperature  and  deposition  rate  were  compared  with  data. 
Good  agreement  was  observed  between  the  model  predictions  and  measured 
data. 

10.  Large  Eddy  Simulations  (LES)  of  the  flow  of  SFg  in  a  rectangular  channel  was 
done  to  resolve  the  inherent  unsteadiness  (if  any)  in  the  flow  of  supercritical 
fluids.  Results  showed  that  unsteady  behavior  (induced  solely  by  gravity)  was 
present  only  for  low  Reynolds  number  flows  and  horizontal  flow 
configurations.  For  higher  Re  and  upward  /downward  flow  configurations, 
the  model  converged  to  a  steady  state  solution. 

11.  Heat  transfer  experiments  were  performed  at  the  University  of  Iowa  using 
supercritical  SF6  as  a  simulant  fluid.  Wall  temperatures  and  exit  fluid 
temperatures  were  measured  for  various  flow  rates  and  furnace  temperatures 
corresponding  to  conditions  ranging  from  the  sub-critical  state  to  the 
supercritical  state. 
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12.  The  SF6  heat  transfer  experiments  were  simulated  using  the  computational 
model.  Good  agreement  was  obtained  over  a  range  of  conditions. 

The  details  of  the  above  tasks  are  discussed  in  this  report. 

1.9  Summary  of  Current  Report 

Section  2  describes  the  models  and  the  governing  equations.  Sections  3  and  4 
describe  the  experimental  work  performed  (at  the  University  of  Iowa  and  Wright 
Laboratory)  during  this  study.  Section  5  presents  the  results  for  the  validation  of  the 
models  for  thermophysical  properties.  The  validity  of  different  turbulence  models 
for  application  to  supercritical  fluid  flow  is  discussed  in  Section  6.  Section  7  presents 
the  results  for  flow  and  heat  transfer  for  supercritical  C02  in  the  laminar  and 
turbulent  regimes.  The  thermal  stability  models  for  aviation  fuels  are  assessed  for 
generality  and  accuracy  and  the  results  are  described  in  Section  8.  The  benchmark 
validation  study  of  the  models  for  wall  deposition  is  described  in  Section  9.  The 
results  from  the  benchmark  validation  study  for  flow  and  heat  transfer  of  SF6  in  a 
tube  are  discussed  in  Section  10.  Section  11  presents  the  results  for  the  LES  study  of 
supercritical  SF6in  a  channel.  Section  12  concludes  with  the  summary  of  the  Phase 
II  work  and  describes  plans  for  future  work. 
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2. 


DESCRIPTION  OF  MODELS 


This  section  describes  the  governing  equations,  the  models  for  thermophysical 
properties  of  supercritical  fluids  flow,  turbulence  and  thermal  stability  of  aviation 
fuels. 


2.1  Governing  Equations 

The  governing  equations  for  the  flow,  heat  transfer  and  chemistry  of  supercritical 
fluids  are  as  follows  [16]: 

Mass  Conservation: 


JtVp»  =  «  (1) 

where  p  is  the  fluid  density  and  u  the  velocity  vector. 

Momentum  Conservation: 

+  V-puu  =  -Vp  +  V-x+pg  (2) 

T  =  peff(Vu  +  VuT-|V-ul) 

where 


M-eff  =  M-  +  fj-t  (3) 

where  p  is  the  hydrodynamic  pressure,  (i  the  dynamic  viscosity,  and  g  the  gravity 
vector.  [it  is  the  turbulent  viscosity  computed  from  the  turbulence  model  and  I  is 
the  unit  tensor. 
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Energy  Conservation: 


+  V  '  Puh  =  V  '  ^eff 


:VT  +  f  +  .: 


Vu 


(4) 


where 


A,g££  —  X  +  At- 


(5) 


where  h  is  the  static  enthalpy,  X  the  thermal  conductivity,  and  x  the  shear  stress 
tensor.  A,t  is  the  turbulent  conductivity  computed  from  (it  using  an  appropriate 
turbulent  Prandtl  number  (Prt). 


Turbulent  Kinetic  Energy  (k): 


^  +  V-pak  =  V 


Vjtt 

<5, 


Vk 


+  p  (P  -  e  -  D) 


Dissipation  Rate  (e): 


^  +  V.pue  =  V 


M-eff 


Ve 


pe^ 


+E 


(6) 


(7) 


The  standard  values  for  the  constants  appearing  in  the  turbulence  model  are: 

Cu  =  0.09 
Ce  =  1.44 
C^=  1.92 

ak  =  1.0 

ae  =  1.3 

The  production  rate  of  turbulence  is  given  by: 

P  =  ^Vu  •  (Vu  +  Vu7)  (8) 
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Species  Conservation: 


dpY; 

“St1  +  V  ‘  PuYi  =  v  •  DiVYi  +  ,  i  =  1, Nsp  (9) 

where  Yj  is  the  species  mass  fraction,  Dj  is  the  mass  diffusivity,  d)j  are  the  reaction 
rate  terms  (computed  as  shown  below)  and  Ngp  is  the  number  of  species. 

For  a  general  chemical  reaction  with  3rd  body  [17]: 


N  Nsp 

X  v..A.  +  M  «•  X  v..A.  +  M  j  =  1, N  . 

i  =  i  i]  i  i  =  i  i)  i  J  st 


(10) 


the  reaction  rate  can  be  expressed  in  the  form: 


=  (|  Cij[A  Jb.  }{KfJ.rf  [A  J><  -  Kb,,^  [AJ»; 


(11) 


where  vjj,  vjj,  Ap  [A  J  are  forward  and  reverse  stoichiometric  coefficients,  molecular 
symbol  and  concentration  of  ith  species,  respectively  and  Nst  is  the  number  of 
reaction  steps.  Qj,  bjj  are  coefficients  of  the  3rd-body  efficiency.  Kf  j,  Kbj,  a  ”,  a--  are  the 
rate  constants  and  concentration  exponents  of  reaction  rate,  and  the  Kf .  ^ 
be  expressed  in  Arrhenius  form: 

Ky^TAexp^)  (12) 


orKb  j),  ran 


For  reactions  obeying  mass  action  law: 


(13) 
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The  reverse  rate  constant  is  derived  from: 


K 


b'i~  K 


eq,j 


(14) 


with  equilibrium  constant: 


(  Po 
UuTj 


RUT 


(15) 


where  pQ  is  1  atm  and  Gibbs  free  energy,  =  ht  -  Tsj. 


The  final  form  of  production  rate  for  ith  species  is  obtained  as: 


(16) 


This  value  is  substituted  into  the  transport  equation  for  species  i. 

2.2  Models  for  Thermophvsical  Properties 

The  properties  such  as  p,  (I,  X,  and  r\  in  the  above  equations  have  to  be  calculated  for 
the  supercritical  regime.  Ideal  gas  approximations  are  not  valid  due  to  steep 
variations  (in  these  properties)  across  the  critical  point.  The  following  sections 
present  a  detailed  formulation  for  the  estimation  of  these  properties. 

Calculation  of  Viscosity:  The  low  pressure  viscosity  is  calculated  from  the 
expression  by  Chung,  et  al.  [18] 


"n 


40.785 


FC(MT)1/2 

V273^, 


(17) 
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T|°  =  low  pressure  viscosity,  |iP 
M  =  Molecular  Weight,  g/mol 
T  =  Temperature,  °K 
Vc  =  critical  volume,  cm3/mol 
Q,  is  a  collision  integral  of  the  form: 

n,  =  [a  (T'H  +  c  [exp  (-DT*)]  +  E  [exp  (-FT’)]  (18) 

where 

T*  =  kT/e,  A  =  1.16145,  B  =  0.14874 
C  =  0.52487,  D  =  0.77320,  E  =  2.16178,  F  =  2.43787 
k  =  Boltzmann  Constant 

e  =  characteristic  energy  of  interaction  between  molecules 

The  collision  integral  is  based  on  a  Lennard-Jones  potential  function  of  the  form, 

¥(r)  =  4e[(f)12-(f)6] 

where  a  is  the  collision  diameter  of  the  molecule. 

The  factor  e/k  can  be  expressed  as: 

£  _  Tc 
k  1.2593 

where  Tc  is  the  critical  temperature.  Therefore, 

T*  =  1.2593  Tr 

where  Tr  =  T/Tc  is  the  reduced  temperature. 

The  factor  Fc  is  calculated  as: 

Fc  =  1  -  0.2576  oo  +  0.059035  p4  +  k  (22) 


(20) 


(21) 
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where  to  is  an  acentric  factor  dependent  on  the  polarity  of  the  molecule,  (ir  is  a 
dimensionless  dipole  moment  and  k  is  a  correction  for  highly  polar  molecules. 

The  high  pressure  viscosity  using  the  Chung  formulation  is  given  by: 


T|=T1 


,36.344  (MTJ172 
Vc2/3 


(23) 


where 

’i* = ^ — (Fc  [(g2)_1 + E6y]} + 1 


y  = 


V  =  molecular  (cm^/  mol) 


_  1  -  0.5y 

1=  (i-y)3 


G2" 


g{[i-exp(~E4y)]/y}  +  G^exp  (E^y)  +  £3 


E^E4  +  E^  +  E^ 


ri“  =  E7  y2  G2  exp  [e8  +  E,  (tT*  +  EI0  (tT* 


The  constants  Ej  through  E10  are  tabulated  in  the  following  table. 
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Table  2-1.  Constants  for  Viscosity  Model 


Ej  =  a;  +  bj  w  +  q  (4  +  d;  k 


/ 

ai 

bi 

ci 

di 

1 

6.234 

50.412 

-51.680 

1189.0 

2 

1.210  x  10'3 

-1.154x10-3 

-6.257x10-3 

0.03728 

3 

5.283 

254.209 

-168.48 

3898.0 

4 

6.623 

38.096 

-8.464 

31.42 

5 

19.745 

7.630 

-14.354 

31.53 

6 

-1.900 

-12.537 

4.985 

-18.15 

7 

24.275 

3.450 

-11.291 

69.35 

8 

0.7972 

1.117 

0.01235 

-4.117 

9 

-0.2382 

0.06770 

-0.8163 

4.025 

10 

0.06863 

0.3479 

0.5926 

-0.727 

4240/1  tl 


The  above  formulation  for  calculating  dense  fluid  viscosities  shows  good  agreement 
with  experimental  data  and  the  errors  are  usually  less  than  5%. 

Calculation  of  Thermal  Conductivity:  The  high  pressure  thermal  conductivity  of  a 
substance  is  given  by  the  expression  of  Chung  et  al.  [18]. 

x=31^T(G_1  +  B^)  +  <iB7y2Ti/2G2  (24) 

X  =  thermal  conductivity,  W/(m  -  °K) 

T)°  =  Low  pressure  gas  viscosity,  N  •  S/m2 
M  =  Molecular  weight.  Kg/mol 
q  =  3.586  x  10~3  M"2  /  V2/3 

The  factor  'P  is  given  as: 

'P  =  1  +  a  {[0.215  +  0.28288a  - 1.061  p  +  0.26665Z]/[0.6366  +  pz  +  1.061ap]}  (25) 
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where 


a  =  X“f '  (cv  =  cP-R) 

P  =  0.7862  -  0.7109oo  +  1.3168  co2 
z  =  2.0  +  10.5  Tr2 

Cv  =  Specific  heat  at  constant  volume 
Cp  =  Specific  heat  at  constant  pressure 
R  =  Universal  gas  constant 


Gi 


1  -  0.5y 

a-y)3 


Bi  r 

^  _TL1_exP("B4y)]  +  B2G1exp(B5y)  +  B3G1 
2  "  Ba  B4  +  B2  +  B3  ~ 

the  coefficients  B:  to  B7  are  tabulated  below: 


Table  2-2.  Constants  for  the  Thermal  Conductivity  Model 

B;  =  aj  +  b  ;co  +  +  dj  k 


1 

2 

3 

4 

5 

6 
7 


2.4166  E+0 
-5.0924  E-1 
6.6107  e+0 
1.4543  E+1 
7.9274  E-1 
-5.8634  E+0 
9.1089  E+1 


7.4824  E-1 
-1 .5094  E+0 
5.6207  E+0 
-8.9139  E+0 
8.2019  E-1 
1.2801  E+1 
1.2811  E+2 


-9.1858  E-1 
-4.9991  E+1 
6.4760  E+1 
-5.6379  E+1 
-6.9369  E-1 
9.5893  E+0 
-5.4217  E+1 


1.2172  E+2 
6.9983  E+1 
2.7039  E+1 
7.4344  E+1 
6.3173  E+1 
6.5529  E+1 
5.2381  E+2 


The  ideal  gas  specific  heat,  Cp/  is  calculated  from  the  JANNAF  data  for  each  species 
as  a  function  of  temperature. 

Calculation  of  Density:  The  Corresponding  States  principle  is  used  to  calculate  the 
density  of  the  fluid. 
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7  _  PV 
RT 


(26) 


where  Z  is  the  compressibility  factor.  Z  is  obtained  as: 

Z  =  z(°)(Tr/Pr)  +  coZW(TrPr)  (27) 

The  function  Z(°)  would  apply  to  spherical  molecules  and  Z(i)  is  a  deviation 
function.  Z(°)  and  ZM  are  tabulated  as  a  function  of  Tr  and  Pr  in  Tables  3-2  and  3-3  of 
Reference  18,  respectively. 

Calculation  of  Enthalpy:  For  a  pure  substance,  the  zero  limit  of  pressure  signifies 
the  enthalpy  of  an  ideal  gas,  h*.  The  real  gas  enthalpy  'IT  at  any  temperature  can  be 
obtained  from  the  following  relationship  as  a  function  of  Tr  and  Pr: 

^=F(0W«  (28) 

where  F(°)  and  F(D  are  tabulated  as  functions  of  Tr  and  Pr  in  Tables  5-2  and  5-3  of 
Reference  18,  respectively.  The  ideal  gas  enthalpy  h*  is  obtained  from  the  JANNAf 
data.  It  is  observed  the  dependence  of  'h'  on  the  temperature  is  highly  non-linear. 
An  iterative  process  is  adopted  to  obtain  the  temperature  from  the  enthalpy  at  any 
given  point  in  the  computational  domain. 

2.3  Models  for  Turbulence 


Most  of  the  turbulence  models  either  explicitly  or  implicitly  assume  the  presence  of 
a  logarithmic  layer  near  the  wall.  Fortunately,  the  logarithmic  layer  is  shown  to 
persist  in  various  wall-attached  flows  such  as  flow  in  adverse  pressure  gradient, 
compressible  flows  etc.  This  is  to  be  expected  because,  far  away,  from  the  wall, 
turbulence  production  dominates  over  dissipation  and  very  close  to  the  wall  (in  the 
laminar  sublayer),  dissipation  is  dominant.  How  ever,  in  some  region  close  to  the 
wall,  the  turbulence  production  will  be  of  the  same  order  as  dissipation.  This  region 
is  called  the  logarithmic  layer.  From  the  law  of  the  wall,  the  turbulent  velocity  and 
temperature  profiles  may  be  expressed  as: 


23 


4241/6 


(29) 


u,  =  tjlnl'Ey  j 

T  (Tw-T)V^ 
+  4„/  pcp 


where 


W 


U  =  —r- 

+  -F 


u 


w 

k  =  0.4034 
E  =  9.7 


/P 


(30) 


p  is  the  fluid  density,  Cp  is  the  specific  heat,  Tw  is  the  wall  temperature,  T  is  the  fluid 

temperature,  qw  is  the  wall  heat  flux,  tw  is  the  wall  shear  stress,  y  is  the  normal 
distance  from  the  wall,  and  u  is  the  fluid  velocity. 

Standard  k-e  Model:  In  the  standard  k-e  turbulence  models,  the  k  and  e  equations 
are  not  integrated  up  to  the  wall.  Instead,  the  shear  stress  and  the  heat  transfer  are 
calculated  using  wall  functions.  Hence  the  center  of  the  first  grid  cell  near  the  wall 
is  required  to  be  approximately  located  in  the  logarithmic  layer  (y+  >30)  and  the 
computational  grid  will  be  relatively  coarse  near  the  wall.  The  wall  functions  used 
in  the  CFD-ACE  code  are  given  below: 


U 


—  C1/4k1/2  =  r-ln 


X  V 
w 


1/2 


(31) 


(T-Tw)CX/<kl/2  P,, 


K 


In 


pEy(cp/2k) 


+  p  -JEZ4  (A) 
rt  sin  7i/4  k  K ; 


1/2^  vpo 

VPr.  APV 


1/4 


(32) 
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where  A  is  the  van  Driest's  constant  take  to  be  26.0.  Pr2  and  Prt  are  the  laminar  and 
turbulent  Prandtl  numbers,  respectively. 

Two-Laver  Turbulence  Model:  As  mentioned  above,  in  the  standard  k-s  model  the 

shear  stress  and  heat  transfer  are  calculated  from  k  and  £  values  at  the  first  grid  point 
near  the  wall  using  wall  functions.  It  has  been  shown  that  these  predictions  are 
sensitive  to  the  location  of  the  near-wall  grid  point  especially  for  supercritical  fluid 
flows.  To  circumvent  these  difficulties,  a  two-layer  model  was  proposed.  In  this 
approach,  the  near-wall  region  is  divided  into  two  layers:  (1)  an  inner  layer  in 
which  the  viscosity  length  scale  is  obtained  form  a  Van  Driest  type  correlation;  and 
(2)  an  outer  layer  in  which  the  standard  k-e  model  is  used.  The  interface  between 
these  two  layers  is  identified  by  a  matching  criteria.  The  major  advantage  of  this 
approach  is  that  substantially  fewer  grid  cells  are  required  in  the  inner  layer  than 
those  required  for  a  low  Reynolds  number  model.  The  mathematical  details  of  the 
two-layer  model  are  described  below  and  the  results  obtained  using  this  model  are 
presented  in  Section  6. 

In  the  outer  layer,  the  standard  k-£  equations  govern  the  transport  of  turbulence  and 
the  eddy  viscosity  is  computed  as: 


M't  ~  e  (33) 

where  is  taken  to  be  0.09.  In  the  inner  layer,  the  molecular  viscosity  is  either 
dominant  or  comparable  to  the  eddy  viscosity.  Hence  the  above  equation  is  replace 
by. 


(34) 

where  is  a  viscosity  length  scale.  It  is  determined  using  the  following  Van  Driest 
type  correlation: 


/n=c/y[1-exp(-¥). 


(35) 
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where  Re  is  a  local  Reynolds  number  based  on  the  turbulent  kinetic  energy  given  by 


(36) 


and  the  Q  is  constant.  To  satisfy  the  log  law  near  the  wall,  Q  is  taken  to  be 


(37) 


where  k  is  the  Von  Karman  constant.  In  the  outer  layer,  the  e  equation  is  replaced 
by  an  algebraic  length  scale  equation: 


.  3/2 

k 


where  the  dissipation  length  scale  is  modeled  as: 


c/y 

e  1+b 
Re 


(38) 


(39) 


In  the  CFD-ACE  code,  the  values  for  constants  a  and  b  appearing  in  the  above 
equations  are  taken  to  be  [19]: 

a  =  50.5  and  b  =  5.3 


The  ratio  between  the  viscous  and  dissipation  length  scales  is  given  by: 

f  _  l\i  _  1-  exp  (-Re/a) 

T»“  "  1/(1 +  5.3 /Re) 


(40) 


Note  that  is  only  a  function  of  Reynolds  number.  The  two-layer  model  is  based 
on  the  premise  that  the  viscosity  length  scale  is  smaller  than  the  dissipation  length 
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scale  inside  the  inner  layer  and  approaches  unity  towards  the  outer  layer.  Though 
several  criteria  have  been  used  to  identify  the  interface  between  the  layers,  the 
location  where  f^  becomes  unity  is  identified  as  the  interface  in  the  CFD-ACE  code. 

Low  Reynolds  Number  Model:  In  the  case  of  near-  or  super-critical  flows,  the 
transport  properties  such  as  molecular  viscosity  and  thermal  conductivity  are  strong 
functions  of  local  conditions.  Consequently,  the  constants  used  in  the  law  of  the 
wall  may  not  be  constants  and  could  vary  for  different  fluids.  Hence  the  standard  k- 

e  model  may  not  be  appropriate  for  supercritical  fluids.  Moreover,  since  the  grid 
will  be  relatively  coarse  near  the  wall,  the  local  property  variation  cannot  be  taken 
into  account.  The  low  Reynolds  number  model  requires  the  grids  to  be  very  fine 
near  the  wall  (y+  <1).  In  the  low  Reynolds  number  model,  the  eddy  viscosity  is 
given  by  [20]: 


k2 

Ft -  pCji  f^  ~ 


where  fM  is  given  by  [20]: 


^  =  1  -  exp  (-0.0115  y+)  (41) 

The  function  fl7  f2,  E  used  in  the  e  equation  are  given  by: 

E  =  -2pp  exp  (-  0.5  y+)  and  fx  =  1.0;  f2  =  1  -  0.22  exp  (~(^f)2)  (42) 


In  the  standard  k-£  model,  the  turbulent  Prandtl  number  is  taken  to  be  constant  (0.7). 
For  supercritical  fluids,  however,  this  assumption  does  not  yield  accurate  results.  It 
has  been  observed  that  the  turbulent  Prandtl  number  close  to  the  wall  increases  due 
to  the  sublayer  "bursting"  phenomenon.  In  order  to  account  for  these  effects,  the 
following  correlation  is  used  in  the  present  calculations  [21]: 


(43) 
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where  Pe  is  the  turbulent  Peclet  number  given  by: 


and  Prhi  is  the  turbulent  Prandtl  number  far  away  from  the  wall,  an  experimentally 
determined  constant  and  C  is  an  experimental  constant.  For  the  results  shown  in 
the  next  section,  C  is  taken  to  be  0.2  and  Pr  is  taken  to  be  0.86. 

2.4  Models  for  Thermal  Stability  of  Fuels 

The  process  of  thermal  decomposition  of  fuels  that  leads  to  surface  deposition 
consists  of  a  series  of  complex  reaction  mechanisms.  The  reaction  initiation 
depends  on  the  temperature  of  the  fuel.  For  fuel  temperatures  less  than  about 
540°K,  the  primary  mechanism  for  deposition  is  the  oxidation  of  fuel  by  the 
dissolved  oxygen,  i.e.,  auto-oxidation.  For  fuel  temperatures  greater  than  750°K, 
thermal  pyrolysis  of  hydrogen  molecules  and  the  scission  of  hydrogen  atoms  from 
carbon  atoms  become  the  dominant  processes.  In  the  intermediate  temperature 
range,  deposition  occurs  through  the  formation  of  precursors  which  decompose  at 
the  wall  to  give  rise  to  the  deposit. 

Thermal  stability  of  jet  fuels  has  received  considerable  attention  over  the  last 
decade.  Recent  studies  have  focused  on  the  identification  and  modeling  of  various 
reaction  mechanisms  for  thermal  decomposition  of  fuels.  Two  models,  the  4-step 
model  and  the  9-step  model  are  described  below. 

4-Step  Model:  Krazinski  et  al.  [22]  proposed  the  following  global  fuel  deposition 
mechanism  for  bulk  fuel  and  surface  reactions  at  the  wall.  The  auto-oxidation  of 
the  fuel  with  the  dissolved  oxygen  can  be  expressed  as: 


ki 

Fuel  +  Oxygen  — »  Precursors  (533°K  <T  <755°K) 

The  precursor  may  either  lead  to  direct  deposition  at  the  wall  or  to  the  formation  of 
solubles  at  higher  temperatures: 
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k2 

Precursors  -»  Solubles(533°K  <T  <755°K) 


The  reactions  at  the  wall  may  be  expressed  as: 


Fuel  +  Oxygen  Deposit  (T  <533°K) 


Precursors  Solubles(533°K  <T  <755°K) 


For  fuel  temperatures  less  than  533°K,  the  first  surface  reaction  controls  the 
deposition  process  without  the  formation  of  precursors.  For  temperatures  between 
533°K  and  755°K,  precursor  formation  (and  transformation  to  solubles)  according  to 
gas  phase  reaction  1  and  2  is  the  controlling  mechanism.  The  deposit  formation  at 
the  wall  occurs  due  to  the  decomposition  of  the  precursor  at  the  wall  according  to 
the  second  surface  reaction.  For  temperatures  greater  than  755 °K,  fuel  pyrolysis 
becomes  the  dominant  process. 

The  rates  for  the  above  reactions  are  written  as: 

k,  =  10n  exp  (-  X 103)  m3 

v  T  1  kg  •  mol 

k2  =  3xl015exp(-iZ4pO!)sec-l 


for  surface  reactions  at  the  wall. 


n  acoj 

D°2  ~W~ 


=  k 


•w 

ks  =  6.34  x  10"5  exp  (-  x  —  )  — — — - 
1  v  T  J  kg  •  mol 


•Sl[Fuel]w[02]w 


m4 


(45) 
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The  precursor  is  assumed  to  decompose  at  the  wall  with  a  reaction  probability  of 
one,  i.e.,  the  precursor  concentration  on  the  wall  can  be  set  to  zero  and  the  precursor 
flux  can  be  directly  translated  into  an  equivalent  deposition  rate,  i.e.. 


Pr]w  =  °  <«) 

In  the  present  study,  transport  equations  were  solved  for  oxygen  and  precursor 
concentration  in  the  field  using  the  rate  expressions  given  above.  Since  oxygen  and 
the  precursors  exist  as  trace  species  in  the  fuel,  their  impact  on  bulk  flow  and  heat 
transfer  characteristics  are  ignored.  Therefore,  the  flow  field  and  the  heat  transfer 
are  calculated  using  the  single  component  (i.e.,  the  fuel)  approximation  and  the 
oxygen  and  precursor  concentrations  are  treated  as  passive  scalars. 

The  total  deposition  rate  is  calculated  as: 


m  = 


-  (Mfuel  +  M02|d02 


atoj 


+  D 


a[Pr] 
Pr_ay“ 


IwJ 


(47) 


where  M  is  the  molecular  weight  and  m  is  the  mass  flux  to  the  surface. 

9-Step  Model:  A  more  advanced  model,  developed  at  Wright  Laboratory  by  Katta  et. 
al.  [23]  was  also  incorporated  into  CFD-ACE.  This  chemistry  mechanism  is 
illustrated  below: 

Reaction  Activation  Temperature  (°K)  Pre-Exponential 

Gas  Phase: 

F  +  02  -4  R00°  16.1xl03  2.5xl010(kmol/m3sec) 

k 

R00°  +  F  4  Solubles  5.04xl03  lxl04(sec_1) 

R00°  +  F  4  P  7.55xl03  8.0xl09(sec_1) 

S 
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R00°  +  F  D 


bulk 


P  +  F  — >  Solubles 


Dbulk  +  F  2Dbulk 


Surface: 


02  +  F  P 


P  D 


wall 


^bulk  ^  ^wall 


5.04xl03 


15.1x10s 

0.0 


6.04x10s 

8.56x10s 

5.04x10s 


2.0xl02(sec'1) 

3.2xl012(sec'1) 

1.0xl0'3(sec'1) 

5.2xl0‘s 

260.0  (m/sec) 

0.80  (m/sec) 


The  deposition  rate  is  calculated  as: 

3D 

^  5n  =  ^07  (PI wall  A  8  exp  (-Ta8/T)  +[Dbulk]wall  A  9  exp  (— Ta9 /x)}  (48) 


During  the  computations,  it  was  found  that  the  above  model  agreed  with  data  only- 
under  certain  restrictive  conditions  (See  explanation  in  Section  8).  In  order  to 

improve  the  generality  of  the  model,  the  surface  reactions  were  recalibrated.  The 
recalibrated  model  is  described  below. 

Recalibrated  9-Step  Model:  The  gas  phase  reactions  are  the  same  as  the  original  9- 
step  model.  The  surface  reactions  are  recalibrated  as  follows: 
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Reaction 

Activation  Temperature  ("10 

Pre-Exponential 

Surface: 

02  +  F  ->  Dwall 

4.03  xlO3 

2.0xl0'5  (m4/kmol-sec) 

^8 

P  Dwall 

8.56xl03 

40  (m/sec) 

k9 

Dbulk  Dwall 

5.04xl03 

0.8  (m/sec) 

and  the  deposition  rate  is  given  by: 


D  -  {K  exP  +  A» exp  ^[P]wan  +  A»  -p  ^  K  J 


wal! 

(49) 


For  JP-5  ,  the  constants  are  recalibrated  as: 


A7  =  6.3  x  10  5  L — 2^- - )  and  A8  =  400  m/ s 

'  Vkmol-secy  0 


These  constants  are  closer  to  the  values  used  by  the  original  Krazinski  model. 
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3. 


EXPERIMENTAL  DEPOSITION  STUDIES  WITH  JET  FUELS 


This  section  describes  the  experiments  (conducted  at  Wright  Laboratory)  with  jet 
fuels  over  a  range  of  operating  conditions.  The  objective  of  the  experiments  is  to 
establish  a  data  base  of  temperature  and  deposition  during  the  heating  of  jet  fuel  to  a 
near  thermodynamic  critical  state. 

To  study  deposition  by  jet  fuel  at  high  temperatures,  fuel  was  pumped  through  a 
long  heated  test  section  at  various  flow  rates,  pressures,  and  heat  fluxes,  and  for 
various  lengths  of  time,  and  the  amount  of  deposition  measured  as  a  function  of 
location  along  the  tube.  Flow  rates  of  100,  200,  and  300  ml/min  were  used,  at 
pressures  of  650  and  400  psig.  (The  ambient  pressure  was  assumed  to  be  standard 
condition  and  the  reduced  pressure  as  typically  higher  than  1.3.)  Heat  flux  was 
controlled  through  the  wall  temperature  of  the  radiant  furnace,  which  was  set  at 
either  2192  or  1832°F  (1200  or  1000°C).  Test  duration  was  typically  5  or  10  hours. 

3.1  Apparatus 

The  experimental  setup  is  shown  in  Figure  3-1.  Air-saturated  Jet  A  2926  was  used  as 
the  fuel  in  most  cases.  The  fuel  was  pumped  from  a  supply  tank  through  a 
pulsation  damper  and  a  0.5  pm  filter  into  the  heated  test  section  of  stainless  steel 
tubing.  The  tube  was  heated  by  a  radiant  heater  having  a  maximum  wall 
temperature  of  2192°F  (1200°C).  After  emerging  from  the  test  section,  the  fuel  passed 
through  a  cooling  heat  exchanger,  a  2  pm  filter  to  collect  particulate  debris  resulting 
from  the  heat  stress,  a  needle  valve  to  control  the  pressure,  and  finally  into  a 
collection  tank.  Pressure  was  measured  at  the  furnace  inlet  and  at  the  back  pressure 
valve  in  psig.  A  three-way  valve  located  at  the  collection  tank  allowed  the  flow  to 
be  diverted  for  flow  rate  measurements. 
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Figure  3-1.  Schematic  of  the  Experimental  Setup 
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Saturation  of  the  fuel  with  either  air  or  nitrogen  was  achieved  by  using  the 
appropriate  gas  from  the  building  supply  and  bubbling  it  through  the  fuel  in  the 
supply  tank.  In  the  case  of  nitrogen,  the  tank  was  also  sealed  except  for  a  small  vent, 
and  left  overnight  prior  to  the  test,  to  ensure  the  removal  of  all  oxygen  from  the 
fuel. 

The  pump  (American  Lewa  Model  EK-01)  used  was  a  reciprocating  pump  with  a 
variable  stroke  length.  The  stroke  length  was  set  to  8-mm  in  all  cases,  and  the  speed 
of  the  motor  varied  with  the  electronic  controller.  The  electronic  package  controlled 
the  motor  (rotational)  speed,  not  the  flow  rate.  As  a  result,  flow  rates  needed  to  be 
measured  manually.  This  was  done  by  diverting  the  exit  flow  from  the  collection 
tank  into  a  graduated  cylinder.  A  manually  pre-determined  calibration  curve  was 
used  to  approximate  the  motor  speed  necessary  to  achieve  a  desired  flow  rate. 

The  chamber  of  the  pulsation  damper  (Liquid  Dynamics  Model  Pulsetwin)  was 
filled  with  nitrogen  to  560  psig  for  the  650  psig  tests,  and  320  psig  for  the  400  psig 
tests.  The  damper  usually  worked  well,  resulting  in  no  noticeable  pressure 
variations  on  the  analog  gauges.  If  the  system  pressure  dropped  below  the  nitrogen 
chamber  pressure,  the  system  pressure  began  to  pulsate  severely.  This  only  occurred 
during  start-up  and  shut-down,  however. 

The  tubing  used  in  all  cases  was  thick-walled  AISI  316  stainless  steel  tubing  having 
1/8"  O.D.,  and  0.055"  I.D.,  although  the  second  heat  exchanger  used  thin  walled 
tubing,  having  an  I.D.  of  0.085".  An  important  note  is  that  the  200  and  300  ml/min 
tests  were  not  evenly  heated,  due  to  multiple  passes  through  the  furnace. 
Continuous  tubing  was  used  in  all  cases,  but  the  U-shaped  bends  in  the  tube  rested 
just  inside  the  collars  at  the  ends  of  the  furnace,  where  the  heating  was  not  as 
intense. 

Test  section  wall  temperature  measurements  were  made  with  K-type 
thermocouples  spot-welded  onto  the  outer  surface  of  the  tube  at  known  distances 
from  the  beginning  of  the  tube.  The  exit  fuel  temperature  was  measured  with  a 
thermocouple  probe  housed  in  a  T-connector  at  the  end  of  the  test  section.  The  final 
wall  temperature  thermocouple  was  spot-welded  onto  the  outer  surface  of  this  T- 
connector.  The  beginning  of  the  tube  was  located  about  two  inches  outside  one  end 
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of  the  furnace,  and  in  all  but  the  first  couple  of  tests,  the  exit  T-connector  was  located 
completely  within  the  collar  at  the  other  end  of  the  furnace.  The  furnace  was  a 
radiant  heater  (Lindberg  Model  #55662)  with  a  maximum  wall  temperature  of 
2192°F  (1200°C).  To  check  the  symmetry  of  heating,  the  wall  temperature  profile  of  a 
tube  carrying  only  a  slow  flow  of  nitrogen  was  measured.  The  results  are  shown  in 
Figure  3-2,  along  with  the  dimensions  of  the  furnace.  The  results  demonstrate  that 
the  heating  is  significantly  lessened  within  about  five  inches  of  the  furnace  collars, 
but  it  is  otherwise  quite  even.  It  should  be  mentioned  that  there  existed  errors  in 
furnace  temperature  readouts  as  indicated  by  Figure  3-2  that  the  measured  tube  wall 
temperature  was  higher  than  the  furnace  (wall)  temperature. 


Distance  along  tube,  in. 

Figure  3-2.  Furnace  Temperature  Profile 
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The  heat  exchanger  used  in  the  first  two  tests  consisted  of  several  feet  of  large  coaxial 
tubing,  fitted  around  the  primary  tubing,  and  through  which  water  was  pumped. 
This  arrangement  was  not  efficient,  and  was  replaced.  The  second  heat  exchanger 
consisted  of  roughly  15  feet  of  tubing  which  was  coiled  to  fit  inside  a  reservoir  of 
water. 

Carbon  deposition  was  measured  by  carbon  bum-off  (Leco  Model  RC-412  Multiphase 
Carbon  Determinator).  A  O.OOlg  sample  weight  and  an  oven  temperature  of  650°C 
were  used.  The  accuracy  of  the  determinator  is  specified  by  the  vendor  as  ±3%  of  the 
answer,  or  ±0.02%  of  the  sample  weight,  whichever  is  greater.  Calibration  samples 
were  occasionally  rim,  and  yielded  agreement  within  a  few  percent. 

3.2  Procedure 


After  installing  the  test  section  and  fresh  filters,  and  attaching  the  thermocouples, 
fuel  was  pumped  through  the  system  at  the  operating  pressure  to  check  for  leaks. 
Once  the  system  was  determined  to  be  leak-free,  the  furnace  was  closed  and  layers  of 
insulation  were  placed  over  the  openings  at  the  ends  of  the  furnace. 

When  the  test  was  ready  to  begin,  the  pump  was  turned  on  and  set  to  the  pre¬ 
determined  speed,  and  the  needle  valve  was  adjusted  until  the  desired  pressure  was 
reached.  The  furnace  was  then  turned  on  and  allowed  to  warm  up  for  10  minutes. 
From  this  point,  the  equipment  was  left  to  operate  for  the  desired  test  time.  During 
this  time,  the  temperatures  and  flow  rate  were  monitored;  the  pressure  system  was 
adjusted  or  the  fuel  supply  tank  replenished,  if  it  became  necessary. 

Roughly  20  minutes  were  needed  for  the  temperatures  to  stabilize  after  the  furnace 
was  turned  on.  In  an  attempt  to  account  for  this  initial  warm-up  time,  10  minutes 
were  added  to  the  duration  of  each  run,  so  that  a  5  hour  run  would  last  5  hours  and 
10  minutes  from  the  time  the  furnace  was  turned  on,  to  the  time  it  was  turned  off. 
Also  of  note  is  that  tests  longer  than  5  hours  were  usually  divided  into  several  runs 
on  successive  days,  with  no  changes  made  to  the  system  except  for  the  replacement 
of  the  exit  filter. 
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After  the  furnace  was  turned  off,  the  insulation  covering  the  openings  was 
removed,  and  the  furnace  was  allowed  to  cool  to  roughly  950°F  (which  took  5  to  10 
minutes),  at  which  time  the  furnace  was  opened,  and  the  temperatures  began  to 
drop  very  rapidly.  When  the  furnace  walls  dropped  to  roughly  350°F  (at  which 
point  the  tube  walls  were  roughly  100°F,  and  the  fuel  temperature  was  even  lower) 
the  fuel  flow  was  stopped,  and  nitrogen  from  the  building  supply  was  blown 

through  the  system  to  purge  it  of  fuel,  and  eliminate  any  further  reactions  in  the 
tube. 

After  the  test  section  was  sufficiently  cool,  it  was  removed,  sanded,  and  cut  into  2- 
inch  sections.  The  sections  were  then  rinsed  in  hexane  to  remove  any  remaining 
fuel,  and  dried  in  a  vacuum  oven  in  preparation  for  being  analyzed  in  the  carbon 
determinator.  The  filters  were  similarly  rinsed  and  dried.  Occasionally,  small  flakes 
of  black  material  would  be  noticed  in  the  exit  filter  housing.  These  were  included, 
to  the  best  effort,  in  the  filter  analysis. 

3.3  Results  and  Discussion 

Surface  deposition  data  are  tabulated  in  Appendix  A.  The  exit  filter  deposition  rates 
are  presented  in  Appendix  B.  Wall  temperature  data  are  tabulated  in  Appendix  C. 
Profiles  of  deposition  and  temperature  along  the  tubes  are  presented  in  graphic  form 
beginning  with  the  100  ml/min  tests  in  Figures  3-3  and  3-4.  The  100  ml/min  tests 
all  show  a  similar  profile.  There  is  an  early  region  of  relatively  light  deposition, 
followed  by  a  large  peak,  and  a  decline  to  very  light  deposition  again.  These  two 
regions  apparently  represent  two  different  modes  of  deposition,  and  the  end  of  the 
deposition  presumably  occurs  where  the  oxygen  has  been  completely  consumed. 
The  exception  to  the  pattern  is  the  low  heat  flux  run  with  the  furnace  set  at  1832°F 
(or  1000°C  noted  in  the  figure)  instead  of  2192°F.  In  this  case,  the  early  region  has 
disappeared,  and  the  deposition  peak  has  moved  downstream,  reflecting  the 
increased  time  required  to  heat  the  fuel  to  the  temperature  at  which  the  heavy 
deposition  begins.  Two  tests  not  shown  ended  with  a  power  failure  to  the  pump. 
The  furnace  was  quickly  opened  and  the  system  purged  with  nitrogen,  but  the  tube 
became  very  hot  nevertheless.  The  deposition  profiles  in  these  cases  were 
surprisingly  similar  in  general  shape  to  the  others,  but  they  are  a  bit  more  jagged.  It 
is  believed  that  the  pulsation  damper  provided  a  limited  amount  of  fuel  flow  even 
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The  test  at  400  psig  (#7/94-10),  resulted  in  a  deposition  profile  nearly  identical  to  the 
corresponding  test  at  650  psig  (#7/94-6),  except  that  the  peak  was  shifted  downstream 
about  8  inches,  e.g.,  see  Figure  3-5.  The  temperature  profile  (Figure  3-6)  shows  that 
the  650  psig  test  reached  higher  temperatures  than  the  400  psig  test.  The  higher 
temperature  may  be  responsible  for  the  earlier  peak  observed  in  the  650  psig  test  case 
(#7/94-6).  If  this  is  the  case,  then  the  results  suggest  that  the  pressure  (or  criticality) 
effect  on  deposition  profile  is  not  significant.  The  nitrogen-saturated  test  (#7/94-11) 
was  conducted  under  identical  conditions  to  the  previous  test,  except  that  the  fuel 
was  sparged  with  nitrogen  rather  than  air.  Under  this  condition,  the  carbon 
deposition  is  orders  of  magnitude  below  that  of  the  air-saturated  test  case  as  shown 
in  Figure  3-5;  therefore,  it  can  be  deduced  that  the  dissolved  oxygen  in  jet  fuel  is  a 
major  contributor  to  carbon  deposition. 

Certain  similarities  and  differences  can  be  seen  between  the  200  ml /min  tests  and 
the  100  ml/ min  tests,  e.g.,  by  comparing  Figure  3-7  to  3-3.  First,  the  deposition  is 
separated  more  clearly  into  two  regions:  the  early  region  of  light  deposition,  and  the 
second  region  composed  of  a  large  peak.  In  this  case  the  second  region  has  begun  to 
split  as  well,  showing  two  peaks,  and  suggesting  different  deposition  mechanisms. 
In  the  two  cases  of  high  heat  flux,  the  second  peak  occurs  shortly  after  the  wall 
temperature  has  passed  1000°F,  (see  Figure  3-8)  and  lies  near  a  wall  temperature 
peak  of  about  1100°F.  Pyrolitic  deposition  is  one  possible  explanation  for  this  second 
peak,  but  it  doesn't  explain  the  second  peak  in  the  low  flux  case.  Also,  the  second 
region  has  moved  a  considerable  distance  downstream.  This  is  presumably  due  to 
the  increased  distance  required  for  the  fuel  to  heat  up  to  the  necessary  temperature, 
although  there  is  no  clear,  single  temperature  at  which  the  deposition  begins.  The 
reason  the  five-hour  deposition  lies  farther  downstream  than  the  10-hour 
deposition  can  be  found  in  the  temperature  profiles.  The  200  ml/min  tests  were 
passed  through  the  furnace  twice,  and  in  the  case  of  the  five-hour  test,  the  bend  in 
the  tube  projected  significantly  into  the  cooler  space  of  the  furnace  collar,  thus 
increasing  the  distance  required  to  heat  the  fuel.  In  both  cases,  the  peak  begins  after 
the  fuel  has  passed  through  the  bend  and  where  the  wall  temperature  rises  through 
roughly  900°F. 
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Figure  3-5.  Deposition  Profiles  at  System  Pressure  400  or  650  psig 
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Figure  3-6.  Wall  Temperature  Profiles  at  System  Pressure  400  or  650  psig 
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Distance  along  tube,  in. 


In  comparison  to  the  slower  tests,  the  first  region  of  deposition  in  the  300  ml /min 
test  has  almost  disappeared,  cf.  Figure  3-9.  Unfortunately,  even  three  passes  through 
the  furnace  was  not  enough  for  this  flow  rate,  and  the  test  section  was  not  long 
enough  to  collect  much  of  the  heavy  deposition.  However,  the  deposition  began 

early  enough  to  show  that  it  began  at  roughly  the  point  where  the  wall  temperature 
exceeded  900°F. 

Summarizing  the  deposition  profiles,  it  seems  that  several  different  mechanisms 
are  at  work,  accounting  for  different  regions  of  deposition.  A  rough  estimate  for  the 
minimum  temperature  required  for  all  deposition  encountered  here  is  900°F.  In 
most  cases,  deposition  occurs  only  where  the  wall  temperature  is  at  or  above 
roughly  900°F.  The  first  region  of  deposition  increases  in  strength  with  increasing 
temperatures  (decreasing  flow  rate),  and  is  spread  out  along  the  tube,  suggesting  that 
this  mechanism  may  consist  of  deposits  forming  at  the  wall,  due  to  wall 
temperature  alone.  The  second  region,  possibly  created  by  more  than  one 
mechanism,  begins  after  the  fuel  has  passed  through  a  certain  thermal  history. 
Perhaps  in  this  case  a  minimum  bulk  fluid  temperature  is  required,  which  would 
suggest  that  the  deposits  form  in  the  fluid  and  are  transferred  to  the  walls  by 
convection  (e.g.,  by  turbulence).  All  deposition  vanishes  with  the  removal  of 
dissolved  oxygen,  as  demonstrated  by  the  de-oxygenated  run. 

Following  the  profile  graphs  are  graphs  showing  the  total  surface  deposition  and 
filter  deposition  as  a  function  of  test  duration.  If  the  two  tests  which  ended  in  a 
power  failure  are  neglected,  as  well  as  the  low  heat  flux  test,  the  100  ml/min  tests 
show  a  linear  surface  deposition  rate  (Figure  3-10),  beginning  after  a  two  hour  lag 
time  which  is  probably  the  time  required  to  coat  the  metal  completely  with  deposits. 
This  indicates  that  deposits  form  on  metal  much  less  readily  than  on  other  deposits, 
and  that  the  nature  of  the  deposits  doesn't  change  in  time.  Since  two  data  points 
can  always  be  used  to  form  a  straight  line,  the  200  ml/min  data  is  inconclusive.  If 
such  a  line  is  drawn,  however,  it  agrees  nicely  with  the  100  ml/min  data.  The 
deposition  in  this  case  begins  after  about  2  1/2  hours,  and  rises  at  roughly  twice  the 
rate  as  in  the  100  ml/min  tests.  The  300  ml/ min  data  point  is  not  included,  since 
the  deposition  profile  shows  that  the  fuel  at  the  exit  had  not  yet  consumed  all  of  the 
dissolved  oxygen.  The  two  tests  indicated  with  an  asterisk  include  estimated  data 
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points  in  the  deposition  profile.  These  data  points  are  missing  because  the  2-inch 
tube  sections  which  would  have  provided  them  were  saved  for  surface  morphology 
and  composition  studies.  Three  sections  were  saved  in  each  case,  and  the  deposition 
on  them  was  estimated  by  taking  the  average  of  the  two  adjacent  sections. 


Distance  along  tube,  in. 

Figure  3-9.  Wall  Temperature  and  Deposition  Profile  at  Flowrate  300  ml/min 


Duration,  hr 

Figure  3-10.  Surface  Deposition  as  a  Function  of  Test  Duration  (Furnace 
Temperature  at  2192°F) 


Exit  filter  deposition  as  a  function  of  test  time  is  shown  in  Figure  3-11.  A  general 
increase  of  deposition  with  time  can  be  seen,  but  no  clear  pattern  can  be  discerned. 
The  exit  filter  deposition  occurred  after  the  fuel  had  passed  through  the  15-ft  length 
of  the  heat  exchanger,  and  is  probably  influenced  by  factors  other  than  the  controlled 
test  conditions.  One  clear  influence,  however,  is  the  amount  of  oxygen 
consumption  in  the  test  section.  In  the  case  of  the  300  ml/ min  test,  which  did  not 
undergo  complete  oxygen  consumption,  a  7  pm  filter  was  used,  rather  than  the 
usual  2  pm  filter,  and  large  amounts  of  deposits  were  collected  anyway.  In  the  case 
of  the  two  lower  flux  tests,  the  exit  filters  collected  so  much  deposition  that  they 
began  to  clog  every  few  hours,  after  which  the  pressure  drop  across  the  filter,  as 
measured  by  the  difference  between  the  two  pressure  gauges,  would  rise  steadily, 
and  fairly  quickly. 

Finally,  the  behavior  of  the  pressures  and  temperatures  should  be  noted.  After  the 
initial  brief  temperature  rise  following  the  activation  of  the  furnace,  all 
temperatures  would  generally  decline  slightly  for  about  1-1/2  to  2  hours.  After  this, 
temperatures  in  the  first  half  of  the  tube  would  stay  fairly  constant,  whereas 
temperatures  farther  on,  where  deposition  was  usually  occurring,  rose  steadily.  A 
typical  example  of  this  is  shown  in  Figure  3-12.  The  distances  given  are  distances 
from  the  beginning  of  the  tube.  Temperatures  would  also  vary  somewhat  from  day 
to  day,  as  well  as  with  varying  pressure.  The  temperatures  shown  in  the  profile 
graphs  are  the  temperatures  at  the  end  of  test.  Temperatures  well  within  the 
furnace  would  also  frequently  exhibit  a  small  oscillation,  with  a  period  of  7-8 
seconds,  and  an  amplitude  of  several  degrees  Fahrenheit.  In  some  of  the  later  tests, 
thermocouples  would  sometimes  give  unreasonably  high  temperatures,  or  none  at 
all.  In  some  cases  this  was  due  to  poor  contact  or  complete  loss  of  contact  between 
the  thermocouple  wires  and  the  test  section.  Also,  the  furnace  heating  profile  graph 
shows  that  test  section  wall  temperature  measurements  were  significantly  higher 
than  the  furnace  wall  temperature.  It  is  not  known  how  the  test  section 
measurements  could  be  in  error,  and  the  only  explanation  available  is  that  the 
furnace  thermocouple  may  have  been  pushed  a  small  distance  into  the  insulation  of 
the  furnace  wall,  and  gave  lower  temperature  readings  as  a  result.  Also,  the  high 
temperatures  in  the  furnace  caused  significant  corrosion  damage  to  the 
thermocouple  wires,  and  this  may  have  contributed  to  the  erroneous 
measurements.  The  pressure  was  generally  quite  stable.  It  became  somewhat 
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difficult  to  control,  however,  in  the  case  of  high  flow  or  low  heat  flux,  which  caused 
very  large  build-up  of  deposition  in  the  exit  filter. 

3.4  Summary 

The  results  of  these  tests  of  supercritical  fuel  deposition  show  several  general 
features.  First,  deposition  depends  on  dissolved  oxygen  and  cannot  occur  without  it. 
Next,  pressure  has  a  relatively  insignificant  effect  on  deposition  over  the  range  of 
test  conditions  examined  in  this  study  (which  have  a  reduced  pressure  higher  than 
1.3).  The  dominant  factor  is  temperature,  both  tube  wall  temperature  and  bulk  fluid 
temperature.  Heat  flux  and  flow  rate  have  an  effect,  but  their  effect  is  not  obvious 
apart  from  their  influence  on  the  temperatures.  The  effect  of  test  time  is 
straightforward.  Deposition  continues  to  build  up  as  time  goes  on,  with  a  very 
linear  rate  in  the  case  of  surface  deposition.  Due  to  the  complexity  of  deposition, 
and  the  large  number  of  factors  involved,  these  results  are  not  easily  generalized, 
and  much  more  work  needs  to  be  done  before  a  comprehensive  understanding  of 
deposition  mechanisms  is  gained. 
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Figure  3-11.  Exit  Filter  Deposition  as  a  Function  of  Test  Duration  (Furnace 
Temperature  at  2192°F) 


80  inches 


60  inches 


1 0  inches 


40  inches 

5  inches 

Time,  min 

Figure  3-12.  Wall  Temperatures  as  a  Function  of  Test  Duration,  Experiment  #7/94-1 


4. 


EXPERIMENTAL  HEAT  TRANSFER  STUDIES  WITH  A  SIMULANT  FLUID 


This  section  describes  the  experiments  performed  at  the  University  of  Iowa  with  a 
simulant  fluid  (SF6). 

4.1  Use  of  a  Simulant 


This  study  was  aimed  at  exploring  the  heat  transfer  characteristics  of  a  supercritical 
fluid,  and  evaluating  whether  or  not  pressure  oscillations  occur  during  the  heating, 
without  the  complication  of  multiple  components,  chemical  reactions,  and 
deposition.  As  a  simulant  for  other  more  complex  and  more  hazardous  fluids,  such 
as  jet  fuel,  the  working  fluid  chosen  was  sulfur  hexaflouride  (SF6).  This  was  chosen 
for  its  relatively  low  critical  properties  (45.54°C  and  37.6  bar),  and  its  lack  of  safety 
hazards. 

4.2  Experimental  Setup 

A  diagram  of  the  experimental  setup  is  shown  in  Figure  4-1.  SF6  is  pumped  from  a 
small  storage  tank  at  20.9  bar  (the  vapor  pressure  of  SF6  at  room  temperature)  by  a 
reciprocating  pump.  There  is  no  direct  measurement  of  the  flow  rate,  but  a 
theoretical  value  can  be  calculated  from  the  speed  of  the  pump  motor,  which  is 
given  by  the  electronic  controller,  and  this  value  can  then  be  modified  by  an 
efficiency  factor  supplied  by  the  vendor.  In  the  earlier  experiments  using  jet  fuel, 
the  measured  flow  rates  matched  the  calculated  values  to  within  about  8%.  The 
experimental  setup  also  includes  a  pulsation  damper. 

After  leaving  the  pulsation  damper,  the  fluid  passes  into  the  test  section,  which  is 
heated  by  a  radiant  furnace.  The  furnace  (Mellen  model  TC-12)  has  a  heating  length 
of  66  cm  (26"),  and  a  maximum  wall  temperature  of  1200°C.  The  wall  temperature 
is  measured  with  a  thermocouple  installed  in  the  interior.  A  diagram  of  the  furnace 
is  shown  in  Figure  4-2.  The  heated  test  section  consists  of  blackened  stainless  steel 
tubing  having  an  outer  diameter  of  6.4  or  9.5  mm  (1/4"  or  3/8"),  and  an  inner 
diameter  of  4.6  or  7.9  mm  (0.18"  or  0.31").  Five  K-type  thermocouples  are  spot- 
welded  along  the  length  of  the  test  section  to  measure  the  tube  surface  temperature. 
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Figure  4-1.  Experimental  Setup 
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Figure  4-2.  Furnace  Dimensions 


A  back  pressure  valve  near  the  furnace  exit  allowed  the  system  pressure  to  be 
regulated.  The  hot  fluid  from  the  furnace  is  then  passed  through  a  water-bath  heat 

exchanger  to  cool  it  down  to  room  temperature,  and  finally  back  into  the  storage 
tank.  ° 


Direct  measurements  were  made  of  the  fluid  pressure  and  temperature,  as  well  as 
tube  temperature  and  furnace  wall  temperature.  The  system  pressure  was  measured 
with  a  static  pressure  transducer  (Omega  model  PX  800-900GV)  at  the  furnace  inlet, 
and  a  second  transducer  (PCB  Piezotronics  model  113A24)  at  the  furnace  exit  which 
measures  dynamic  behavior.  Bulk  fluid  temperature  is  measured  with  fine-wire 
thermocouples  inserted  into  the  flow  at  the  furnace  inlet  and  exit.  All 
measurements  were  taken  by  a  computer  based  data  acquisition  system  (Apple  Ilci 
with  National  Instruments  MIO  boards  and  software). 

4.3  Procedure  and  Test  Conditions 

Tests  were  conducted  by  setting  the  flow  rate  and  system  pressure,  and  raising  the 
furnace  temperature  in  increments.  Measurements  of  pressure  and  temperature 
were  taken  continuously  by  the  computer,  and  once  it  had  been  determined  that  the 
system  had  reached  steady-state,  one  set  of  measurements  was  stored  in  the  data  file. 
Once  it  was  discovered  that  this  process  could  be  significantly  accelerated  by  taking 
measurements  of  all  pressures  for  a  given  furnace  setting,  the  tests  were  conducted 
accordingly. 

Two  flow  rates  were  used:  150  ml/min  with  the  7.9  mm  I.D.  tubing,  and  450 
ml /min  with  the  4.6  mm  I.D.  tubing.  These  flow  rates  correspond  to  inlet  Reynolds 
numbers  of  about  2,000  and  10,000,  respectively.  Measurements  were  usually  taken 
at  pressures  from  27.6  to  55.1  bar  (400  to  800  psig),  in  increments  of  6.9  bar  (100  psig), 
in  addition  to  a  measurement  at  the  critical  pressure,  37.6  bar  (531  psig).  The  highest 
furnace  temperature  used  was  1000°C,  which  resulted  in  an  average  heat  flux  of 
about  90  kW/m2  at  a  flow  rate  of  450  ml/min.  Most  tests  were  ended  at  a  furnace 
temperature  of  around  800  C.  Tests  were  conducted  with  both  horizontal  flow  and 
upward  flow. 
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4.4  Experimental  Results 


The  heat  transfer  results  for  the  present  experiments  are  presented  as  function  of 
furnace  wall  temperature,  rather  than  heat  flux.  The  reason  for  this  is  as  follows. 
The  heat  flux  was  controlled  in  these  experiments  through  the  furnace  wall 
temperature.  The  magnitude  of  the  heat  flux  was  not  known,  however.  Rough 
values  could  be  calculated  with  radiation  and  convection  correlations  for  concentric 
cylinders,  but  it  was  found  that  the  heat  flux  was  not  uniform,  due  to  a  number  of 
factors.  Also,  as  with  the  jet  fuel  experiments,  it  is  now  believed  that  the  simple 
model  of  a  smooth,  isothermal  surface  is  inadequate  for  the  real  furnace  wall,  as  the 
wall  is  composed  of  alternating  rows  of  heating  coils  and  ceramic  wall  material. 

In  order  to  compare  heat  transfer  rates  for  various  conditions,  the  rates  were 
calculated  based  on  the  enthalpy  increase  of  the  fluid,  and  plotted  versus  furnace 
temperature.  Figure  4-3a  shows  the  results  for  an  inlet  Reynolds  number  of  10,000, 
various  pressures,  and  the  two  flow  directions.  The  data  has  a  significant  amount  of 
scatter  and  shows  no  clear  dependence  on  pressure  or  flow  direction.  Figure  4-3b 
shows  the  results  for  an  inlet  Reynolds  number  of  2,000.  The  data  for  the  two 
different  flow  directions  were  similar,  and  were  combined  together.  There  is  a  slight 
dependence  on  pressure,  however.  A  comparison  between  the  two  graphs  show  a 
slight,  if  any,  influence  of  Reynolds  number  on  the  heat  transfer  rate.  These  results 
show  that  the  furnace  temperature  is  the  dominating  factor  in  determining  heat 
transfer,  but  other  factors,  such  as  pressure  and  Reynolds  number,  could  influence 
the  rate  for  a  given  furnace  temperature. 
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Figure  4-3.  Heat  Transfer  Rates 


4.4.1  Fluid  Exit  Temperature 

The  enthalpy  increase  of  a  heated  fluid  in  a  simple  forced  flow  situation  is  equal  to 
the  heat  added  to  the  fluid.  Thus,  assuming  a  known  and  constant  pressure,  the 
temperature  rise  of  the  fluid  will  be  directly  related  to  the  heat  transfer  through  the 
enthalpy.  For  a  given  tube  diameter  and  length,  and  a  constant  inlet  temperature, 
the  fluid  exit  temperature  will  be  function  only  of  the  heat  flux.  The  form  of  this 
function  will  be  the  same  as  that  for  fluid  temperature  versus  enthalpy. 

Figure  4-4  shows  the  fluid  exit  temperature  versus  furnace  temperature  for  the 
conditions  of  horizontal  flow  at  the  rate  of  450  ml/min  through  a  4.5  mm  inner 
diameter  tube.  At  very  low  heat  fluxes,  exit  temperatures  at  all  pressures  increase 
steadily  and  there  is  little  difference  between  pressures.  At  intermediate  heat  fluxes 
and  low  pressures,  the  exit  temperature  remains  approximately  constant  and  equal 
to  the  boiling  point  at  each  pressure.  At  supercritical  pressures,  the  exit  temperature 
increases  steadily  with  increasing  heat  flux.  At  high  heat  fluxes,  the  exit 
temperature  increases  rapidly  with  increasing  heat  flux.  Thus  there  is  good 
qualitative  agreement  between  the  observed  exit  temperature  and  the  calculated 
enthalpy  of  the  fluid.  Except  for  the  critical  point  and  the  boiling  points,  which 
agree  adequately,  the  data  cannot  be  compared  quantitatively  because  the  heat  flux  is 
not  known. 

This  relationship  between  fluid  temperature  and  heat  flux  has  led  some  researchers 
of  supercritical  heat  transfer  to  present  their  results  in  the  form  of  graphs  of  wall 
temperature  vs.  bulk  enthalpy.  This  has  the  advantage  of  reducing  two  variables 
(heat  flux  and  distance  along  the  tube)  into  one. 
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Figure  4-4.  Bulk  Fluid  Temperature  at  Exit,  Horizontal  Flow,  Re  =  10,000 
4.4.2  Tube  Wall  Temperature 

The  principal  measurement  of  interest  in  determining  the  effectiveness  of  heat 
transfer  under  supercritical  pressures  is  the  tube  surface  temperature. 
Unfortunately,  it  is  believed  that  local  heat-flux-dependent  errors  were  introduced 
into  the  tube  surface  temperature  measurement  by  the  thermocouples,  especially 
with  the  smallest  tubing.  The  wire  thickness  was  fairly  large  (20  gauge),  and  the 
insulation  covering  the  wires  had  a  tendency  to  become  frayed  a  the  ends,  exposing 
the  wires.  It  is  believed  that  these  factors  led  to  significant  amounts  of  heat 
conduction  to  the  tube  through  the  thermocouple  wires.  These  errors  were 
minimized  as  much  as  possible,  but  were  to  a  certain  extent  unavoidable, 
particularly  with  the  smaller  tubing. 

These  temperature  errors  were  evident  in  jagged  tube  temperature  profiles, 
rendering  them  of  little  value.  However,  if  the  tube  temperature  at  a  given  location 
is  plotted  versus  heat  flux  (bulk  enthalpy),  this  may  substitute  for  a  plot  of  surface 
temperature  versus  distance.  The  closest  surface  thermocouple  was  19  cm  from  the 
exit  fluid  thermocouple,  but  if  these  two  measurements  are  compared  under 
various  conditions,  a  few  qualitative  observations  can  be  made. 
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In  comparing  the  results  from  horizontal  and  vertical  flow,  it  will  be  useful  to  have 
an  estimate  of  whether  buoyancy  effects  are  negligible.  Hall  and  Jackson  [6]  cite  the 
following  experimentally  determined  criterion  for  negligible  buoyancy  effects  in 
supercritical  C02: 

Grb/Reb7<2.4xl0~5 
where  Reb  =  VD/u, 


and  Grb  =  (pb  -  pw)D3g/pb-o2 

Calculating  this  property  ratio  is  difficult  because  of  the  limited  property  data 
available  for  SF6/  but  an  estimate  using  a  flow  rate  of  450  ml/min  and  a  tube  inner 
diameter  of  4.6  mm,  and  property  estimates  at  temperatures  near  the  beginning  of 
the  tube,  yields  a  value  of  about  2  x  10-4,  which  indicates  that  buoyancy  may  play  a 
role.  On  the  other  hand,  a  similar  calculation  using  a  flow  rate  of  150  ml/min,  and  a 
tube  inner  diameter  of  7.9  mm  yields  a  value  of  about  9  x  10'2,  indication  with  much 
more  certainty  that  buoyancy  should  have  an  effect  under  these  conditions.  These 
two  calculations  correspond  to  inlet  Reynolds  numbers  of  10,000  and  2,000, 
respectively. 

Figure  4-5  shows  a  graph  of  the  final  tube  surface  measurement  for  the  two  sets  with 
an  inlet  Reynolds  number  of  10,000.  The  horizontal  flow  graph  corresponds  to 
Figure  4-4.  All  of  the  temperatures  increase  fairly  steadily  for  the  most  part.  An 
important  point  is  that  for  given  flow  conditions,  the  tube  temperature  will 
generally  be  higher  for  higher  pressures.  In  cases  where  the  wall  temperature  and 
fluid  temperature  span  the  pseudo-critical  temperature,  the  increase  of  tube 
temperature  with  pressure  could  be  due  to  the  decrease  of  the  peak  in  specific  heat 
with  increasing  pressure.  Under  other  conditions,  the  general  increase  in  thermal 
conductivity  may  also  play  a  role,  allowing  for  greater  heat  transport  across  the 
boundary  layer.  However  the  viscosity  increases  in  a  similar  manner  to 
conductivity,  and  may  decrease  the  turbulent  transport. 
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Furnace  Temperature  (°C) 
(a)  Horizontal  Flow 


(b)  Upward  Flow 

Figure  4-5.  Tube  Temperature  near  the  Exit,  Re  ~  10,000 
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The  exception  to  the  trend  is  post-crisis  subcritical  heat  transfer.  At  subcritical 
pressures,  the  tube  temperatures  exhibit  a  sharp  increase,  indicating  some  form  of 
heat  transfer  crisis.  Heat  transfer  after  the  subcritical  crisis  is  probably  determined  by 
film  boiling.  If  this  is  the  case,  it  implies  that  the  vapor  layer  against  the  wall 

provides  lower  levels  of  heat  transport  than  the  corresponding  supercritical 
boundary  layers. 

At  high  heat  fluxes  (furnace  temperature  above  about  600°C),  Figure  4-4  shows  that 
at  all  pressures  for  which  data  exist,  the  fluid  exits  the  tube  in  a  gas  or  gas-like  state. 
Figure  4-5a  for  the  same  range  shows  that  the  tube  wall  temperature  generally 
increases  with  increasing  pressure,  in  conjunction  with  the  increasing  fluid 
temperature,  which  is  consistent  with  single-phase  convection.  The  data  for  Pr  = 
0.75,  however,  are  higher  than  all  of  the  rest  in  this  range.  As  indicated  above,  this 
may  be  due  to  poorer  heat  transfer  properties  in  subcritical  vapor-like  fluid.  The 
data  for  Pr  =  0.94  show  a  similar  behavior  at  lower  heat  fluxes,  but  data  at  both 
subcritical  pressures  seem  to  increase  more  slowly  than  the  supercritical  data. 
Perhaps  this  is  due  to  slower  property  variations  with  temperature  in  vapor  than  in 
vapor-like  fluid. 

At  subcritical  pressures,  the  tube  surface  starts  out  at  lower  temperatures  than  it  does 
at  higher  pressures,  but  eventually  jumps  up  to  higher  temperatures,  giving 
evidence  of  a  heat  transfer  crisis  (cf.  Figure  4-5).  The  tube  surface  in  these  low 
temperature  regions  is  somewhat  higher  than  the  boiling  point  and  tends  to 
increase  with  increasing  heat  flux.  This  is  especially  noticeable  at  Pr  =  0.75.  Also,  the 
jump  at  Pr  =  0.94  occurs  at  almost  the  same  time  that  boiling  begins,  suggesting  that 
a  heat  transfer  crisis  occurs  at  low  qualities.  However,  if  nucleate  boiling  is 
occurring  in  the  low  temperature  regions  before  the  heat  transfer  crisis,  the  tube 
temperature  should  remain  approximately  constant  at  the  boiling  point.  Perhaps 

this  discrepancy  is  due  to  increasing  thermocouple  conduction  error  with  increasing 
heat  flux. 

Figure  4-5b  shows  the  results  of  a  test  under  similar  conditions  as  those  for  Figure  4- 
5a,  except  that  the  flow  direction  is  vertically  upward.  Due  to  the  high  degree  of 
scatter  in  the  data,  the  results  are  difficult  to  interpret,  but  they  generally  are  similar 
to  the  horizontal  test.  The  temperatures  for  the  upward  test  are  slightly  higher. 
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however.  This  agrees  with  the  majority  of  existing  data  which  indicate  that  heat 
transfer  is  degraded  in  buoyancy  affected  upward  flow,  but  the  results  in  this  case 
may  also  be  due  to  the  hot  air  in  the  furnace  tending  to  accumulate  at  the  upward 
(exit)  end  of  the  furnace,  where  the  last  surface  thermocouple  was  located. 

Figure  4-6  shows  the  results  of  tests  with  larger  tubing  and  a  slower  flow  rate 
(Re=2,000).  The  results  are  similar  in  general  form,  but  the  temperatures  are  much 
higher  for  a  given  furnace  setting,  due  to  the  lower  flow  rate.  Also,  the  heat  transfer 
crisis  seems  to  be  much  more  gentle.  It  is  not  known  why  this  would  be,  but 
perhaps  it  has  something  to  do  with  buoyancy  forces  causing  bubbles  against  the 
wall  to  rise  and  coalesce,  allowing  a  smoother  transition  from  nucleate  to  film 
boiling.  The  temperatures  in  the  vertical  test  are  noticeably  higher  than  those  for 
the  horizontal  test,  in  agreement  with  the  prediction  of  significant  buoyancy  effects. 
Again,  however,  part  of  this  increase  may  be  due  to  increased  air  temperatures  at  the 
top  of  the  furnace. 

4.4.3  Pressure  Oscillations 

As  with  the  jet  fuel  tests,  no  pressure  oscillations  other  than  those  caused  by  the 
pump  were  observed.  A  calculation  based  on  the  enthalpy  data  for  SF6  indicates  that 
the  highest  average  heat  flux  achieved  in  these  tests  was  about  90  kW/m2.  Virtually 
all  pressure  oscillations  recorded  in  the  literature  with  non-cryogenic  fluids 
occurred  at  heat  fluxes  greater  than  1  MW/m2.  Thus  it  is  likely  that  the  heat  fluxes 
in  the  present  tests  were  insufficient  to  cause  the  oscillations. 
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5.  VALIDATION  OF  MODEL  FOR  THERMOPHYSICAL  PROPERTIES 

The  property  model  for  calculating  the  viscosity,  thermal  conductivity, 
compressibility  factor  (for  density  calculation)  and  the  real  enthalpy  of  the  fluid  were 
incorporated  into  the  code.  This  facilitates  the  computation  of  these  quantities  (as  a 
function  of  reduced  pressure  and  temperature)  for  a  range  of  fluids  in  a  general 
purpose  manner.  Also,  a  database  consisting  of  a  number  of  hydrocarbon  species 
(with  the  respective  properties)  has  been  implemented.  The  code  was  used  to 
perform  calculations  for  a  series  of  problems.  The  results  are  compared  with 
experimental  data  from  the  literature.  The  description  of  each  of  these  cases  is 
provided  below. 

5.1  Properties  of  Carbon  Dioxide 

The  critical  point  properties  of  carbon  dioxide  are: 

pc=  7.38  MPa 
Tc  =  304. 1°K 

The  comparison  for  carbon  dioxide  is  shown  in  Figure  5-1.  The  properties  (density, 
viscosity,  conductivity  and  enthalpy)  of  carbon  dioxide  (predicted  by  the  model)  are 
plotted  as  functions  of  temperature  and  pressure.  It  is  observed  that  the  properties 
predicted  by  the  model  are  in  close  agreement  with  the  data  [7, 18,  24]. 
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5.2 


Properties  of  Propane 


The  critical  point  properties  of  propane  are: 

pc  =  4.25  MPa 
Tc  =  369.8°K 

The  properties  (density,  viscosity,  conductivity  and  enthalpy)  of  propane  (predicted 
by  the  model)  is  plotted  as  a  function  of  temperature  and  pressure  and  is  compared 
with  data  [18,  24]  on  propane  (see  Figure  5-2).  It  is  observed  that  good  agreement 
with  data  is  obtained  over  the  entire  range  of  temperature  and  pressure. 


Figure  5-2a.  Variation  of  Density  with  Temperature  for  Propane 


63 


4241/6 


Figure  5-2d.  Variation  of  Enthalpy  with  Temperature  for  Propane 


5.3  Properties  of  TP-5 

The  critical  point  properties  of  JP-5  are: 

pc  =  22.8  MPa 
Tc  =  685 °K 

The  Comparisons  for  JP-5  are  shown  in  Figure  5-3.  Again,  it  is  observed  that  the 
model  prediction  are  in  reasonable  agreement  with  experimental  data  [25]. 
Deviation  of  model  prediction  for  thermal  Conductivity  from  data  is  observed  for 
the  high  temperature  regime. 
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Viscosity  (Pa-s)  Density  (Kg/n,3) 


JP-5  Properties  at  430  Psia 


SF  Properties 


Temperature  (K) 

Figure  5-4b.  Variation  of  Viscosity  with  Temperature  for  SF, 
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Figure  5-4c.  Variation  of  Thermal  Conductivity  with  Temperature  for  SI 
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Figure  5-4d.  Variation  of  Enthalpy  with  Temperature  for  SF6 


5.5  Summary 

It  is  observed  that  the  model  is  able  to  predict  the  quantitative  trend  in 
thermophysical  properties  as  a  function  of  temperature  and  pressure  for  a  range  of 
fluids.  It  is  especially  important  to  resolve  the  density  variation  since  it  influences 
the  fluid  dynamics  significantly.  Conductivity  variations  in  the  laminar  sublayer 
have  a  significant  influence  on  wall  heat  transfer. 
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6.  ASSESSMENT  OF  TURBULENCE  MODELS 


The  two  layer  model  discussed  in  Section  2  is  expected  to  provide  better  accuracy 
than  the  standard  k-e  model  without  having  to  incur  the  computational  expense  of 
the  low  Reynolds  number  model.  A  series  of  test  cases  were  simulated  to  assess  the 

accuracy  of  the  two-layer  model  against  that  of  the  standard  k-e  model. 
Comparisons  with  experimental  data  are  also  presented.  The  first  set  of  cases 
present  the  comparisons  between  the  model  for  various  fluid  flow  problems.  The 
second  set  of  cases  focus  on  heat  transfer  predictions.  For  all  the  test  cases 
considered,  the  near-wall  grid  clustering  was  maintained  so  that  the  first  grid  point 
away  for  the  wall  lies  in  the  logarithmic  layer  for  the  standard  k-e  model  and  in  the 
laminar  sublayer  for  the  two-layer  and  low  Reynolds  number  models. 

6.1  Comparison  of  Fluid  Flow  Predictions 

Fully  Developed  Channel  Flow:  The  first  case  considered  is  the  fully  developed 
flow  in  a  channel.  The  experimental  data  obtained  by  Laufer  [27]  is  selected  since  it 
provides  mean  velocity  as  well  as  the  turbulence  quantities  at  three  Reynolds 
number:  12300,  30800,  and  61600.  The  Reynolds  number  is  based  on  the  centerline 
velocity  and  half-channel  height.  Only  the  results  for  the  highest  Reynolds  number 
are  presented  here. 

A  schematic  of  the  flow  and  the  computational  domain  is  shown  in  Figure  6-1.  The 
velocity  data  is  integrated  to  obtain  the  average  velocity  which  is  specified  as  the 
uniform  velocity  at  the  channel  inlet.  Thus  the  mass  flow  rates  in  the  experiment 
and  computation  are  maintained  same.  The  channel  is  long  enough  for  the  flow  to 
become  fully  developed.  Though  a  50X40  grid  is  used  for  both  the  models,  the  grid 
is  uniform  in  both  directions  for  the  k-e  model  and  clustered  near  the  wall  for  the 
two-layer  model.  Figures  6-2  and  6-3  show  the  mean  velocity  and  turbulent  kinetic 
energy  across  the  channel.  As  can  be  seen,  the  agreement  between  the  predictions 
and  data  is  very  good  for  the  mean  velocity.  The  peak  turbulent  kinetic  energy, 
however,  is  significantly  overpredicted  by  the  two-layer  model.  This  is  consistent 

with  predictions  of  most  low  Reynolds  number  models  that  integrate  the  k  and  £ 
equations  up  to  the  wall. 
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Figure  6-1.  A  Schematic  of  the  Computational  Domain 
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Figure  6-2.  Variation  of  Mean  Velocity  Across  the  Channel 
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Figure  6-3.  Variation  of  Mean  Turbulent  Kinetic  Energy  Across  the  Channel 

Fully  Developed  Pipe  Flow:  The  experimental  data  obtained  by  Laufer  [28]  is  again 
used  for  validating  the  predictions.  The  grids  used  for  this  case  are  the  same  as 
those  used  in  the  previous  case.  Figures  6-4  and  6-5  compare  the  mean  velocity  and 
turbulent  kinetic  energy  predictions  with  the  data.  Both  the  models  predict  the 
mean  velocity  accurately.  The  peak  turbulent  kinetic  energy  predicted  by  the  k-e 
model  is  lower  than  that  predicted  by  the  two-layer  model.  However,  both  models 
underpredict  the  kinetic  energy  compared  to  experimental  data. 
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Figure  6-4.  Mean  Turbulent  Velocity  as  Function  of  Pipe  Radius 


Figure  6-5.  Mean  Turbulent  Kinetic  Energy  as  a  Function  of  Pipe  Radius 
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Flow  in  Backward-Facing  Step:  The  flow  in  backward-facing  step  has  become  a 
standard  benchmark  problem  to  validate  the  turbulence  models  since  the  1980-81 
AFOSR-HTTM  Stanford  Conference  on  Complex  Turbulence  Flows.  The  data 
obtained  by  Driver  and  Seegmiller  [29]  is  used  for  validation.  A  schematic  of  the 
flow  is  shown  in  Figure  6-6.  The  expansion  ratio  (ratio  of  downstream  channel 
height  to  upstream  channel  height)  is  1.125  and  the  Reynolds  number  based  on  the 
step  height  is  about  36000. 


Figure  6-6.  A  Schematic  of  the  Backward-Facing  Step  Flow 

A  80X60  grid  with  clustering  around  the  step  is  used.  The  clustering  near  the  walls 
is  finer  for  the  two-layer  model.  The  comparison  of  predicted  reattachment  length 
with  the  data  is  shown  below  in  Table  6-1.  Clearly,  the  two-layer  model  prediction 
of  the  attachment  length  is  closer  to  the  data. 

Table  6-1.  Comparison  of  Predicted  Reattachment  Length  Using  Various 


Turbulence  Models  with  Data. 


Model 

Reattachment  Length 

Experimental  Data 

6.0 

Standard  k-e 

5.2 

Two-Layer 

5.6 
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Figures  6-7  and  6-8  show  the  static  pressure  and  skin  friction  along  the  step  wall. 
The  static  pressures  predicted  by  both  the  models  are  identical.  The  discrepancy 
between  the  predictions  and  the  data  may  be  attributed  to  the  underprediction  of 
reattachment  length.  The  predictions  of  skin  friction  by  the  models  differ 
significantly.  The  two-layer  model  predicts  the  minimum  value  of  skin  friction  in 
the  recirculation  region  more  accurately.  However,  it  overshoots  the  data 
downstream  of  the  reattachment  point. 

Flow  in  180  Turnaround  Duct:  Flow  in  turnaround  duct  (TAD)  is  a  popular 
validation  case  for  CFD  codes  with  body-fitted  coordinate  (BFC)  capabilities.  Flow  in 
TAD  is  very  useful  for  validating  turbulence  models  because  the  flow  contains 
boundary  layers  subjected  to  both  convex  and  concave  curvatures.  Thus  both  the 
stabilizing  and  destabilizing  effects  of  curvature  on  turbulence  can  be  studied. 

A  schematic  of  the  TAD  is  shown  in  Figure  6-9.  The  Reynolds  number  based  on 
duct  height  and  the  bulk  velocity  is  about  86000.  The  data  obtained  by  Sandbom  [30] 
is  used  for  validation.  Uniform  flow  is  prescribed  at  the  inlet  and  a  constant 
pressure  is  specified  across  the  exit.  A  80X60  grid  with  clustering  around  the  bend  is 
used.  Figure  6-10  shows  the  static  pressure  variation  along  the  inner  and  outer 
walls.  Pressure  recovery  downstream  of  the  bend  is  captured  better  by  the  two-layer 
model.  Figures  6-11  through  6-13  show  the  mean  velocity  profiles  at  0,  90,  and  180 
around  the  bend.  Note  that  the  two-layer  model  predicts  separation  at  the  inner 

wall  at  180  position,  a  feature  not  captured  by  the  standard  k-e  model. 
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Figure  6-7.  Static  Pressure  Variation  Along  the  Step  Wall 


Figure  6-13.  Mean  Turbulent  Velocity  Profile  at  180 


Flow  in  Rotating  Disk  Cavity:  This  case  tests  the  turbulence  models  for  swirling  and 
rotating  flows.  The  flow  domain,  shown  in  Figure  6-14  consists  of  a  cavity  enclosed 
by  a  stator  and  a  rotor  in  the  axial  direction  and  by  the  shaft  and  a  shroud  in  the 
radial  direction.  Daily  and  Nece  [31]  obtained  mean  velocity  data  for  a  disk  spacing 
ratio  (ratio  of  cavity  width  to  outer  radius)  of  0.0637  at  a  rotational  Reynolds  number 
of  4.4X106.  A  40X60  grid  clustered  near  the  walls  is  employed. 

The  torque  coefficient  is  a  measure  of  the  accuracy  of  the  skin  friction  predictions. 
Table  6-2  below  compares  the  predictions  with  the  data.  Figures  6-15  and  6-16  show 
the  axial  variation  of  radial  and  tangential  velocity  components  at  a  radial  location 
(non-dimensional)  of  0.765.  The  predictions  are  in  good  agreement  with  the  data. 
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Table  6-2.  Comparison  of  Torque  Coefficient  Prediction  With  Data 


Figure  6-15.  Axial  Variation  of  Radial  Velocity 


6.2  Comparison  of  Heat  Transfer  Predictions 


A  point  to  be  noted  is  that  the  analog  of  the  two-layer  model  for  heat  transfer  does 
not  exist.  Prior  applications  and  validations  of  this  model  have  been  limited  to 
isothermal  turbulent  flows  or  to  flows  with  small  temperature  and  property 
gradients  near  the  wall.  Due  to  this  limitation,  the  thermal  wall  functions  of  the 

standard  k-e  model  were  used  (in  this  study)  in  conjunction  with  the  two-layer 
model  to  calculate  heat  transfer  at  the  wall. 

As  discussed  before,  the  two-layer  model  is  expected  to  provide  the  same  level  of 
accuracy  as  the  low  Reynolds  number  model  but  at  a  lower  computational  cost.  In 

other  words,  it  is  expected  to  be  more  accurate  than  the  standard  k-e  model  without 
additional  computational  expense.  Therefore,  by  definition,  one  would  expect  the 
predictions  of  the  two-layer  model  to  lie  in  between  the  prediction  of  the  standard  k- 
e  model  and  the  low  Reynolds  number  k-e  model.  The  results  in  the  previous 
section  show  that  this  is  indeed  true  for  the  prediction  of  pressure  coefficient,  skin 
friction  coefficient  and  mean  velocity  profiles  in  turbulent  incompressible, 
isothermal  flows.  The  comparisons  for  heat  transfer  are  shown  below. 

The  flow  of  CO2  in  a  two  dimensional  channel  is  considered.  The  channel  pressure 
is  7.58  MPa  which  is  above  the  critical  pressure  (Pc=  7.38  MPa)  of  C02.  The  inlet 
temperature  is  288°K  which  is  below  the  critical  temperature  (Tc  =  306.1°K)  of  C02. 
The  channel  wall  is  heated  and  the  fluid  undergoes  transition  to  the  fully 
supercritical  regime  in  the  channel.  The  flow  of  C02  is  upward,  i.e.,  against  the 

direction  of  gravity.  This  problem  is  discussed  in  greater  detail  in  the  next  section 
(see  Figure  7-1). 

I 

Constant  Wall  Temperature  Case:  Figure  6-17a  shows  the  heat  flux  values  predicted 
by  the  three  turbulence  models  for  a  constant  wall  temperature  of  320°K.  the 
computational  grid  for  the  two-layer  model  was  the  same  as  that  for  the  standard  k-£ 
model.  It  is  observed  that  the  heat  transfer  rates  predicted  by  the  two-layer  model  do 
not  lie  in  between  the  values  predicted  by  the  standard  and  low  Reynolds  number  k- 
e  models. 
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Constant  Heat  Flux  Case:  Another  case  was  simulated  with  a  constant  heat  flux  of 
30  kW/m2  applied  to  the  sidewalls.  Again,  it  is  observed  (see  Figure  6-17b)  the 

predictions  of  the  two-layer  model  do  not  lie  between  those  of  the  standard  k-e 
model  and  the  low  Reynolds  number  model.  Two  different  runs  were  performed 
with  the  two-layer  model,  one  with  the  same  grid  resolution  as  the  standard  k-e 
model  and  the  other  with  the  same  grid  resolution  as  the  low  Reynolds  number 
model.  Both  simulations  yielded  results  that  were  less  accurate  than  the  standard  k- 
e  model  predictions. 

The  main  reason  for  this  discrepancy  appears  to  be  the  lack  of  an  appropriate  analog 
for  wall  heat  transfer  in  the  two-layer  model.  This  feature  seriously  limits  the 
application  of  the  two-layer  model  to  turbulent  flow  involving  heat  transfer. 
Therefore  all  turbulent  flow  calculations  presented  in  the  following  sections  of  this 
report  were  done  using  the  low  Reynolds  number  model.  Although  the 
computational  expense  of  this  model  is  relatively  high,  its  ability  to  resolve  near 
wall  property  variations  as  well  as  the  lack  of  any  limiting  assumptions  (such  as 

wall  functions)  makes  it  a  good  candidate  for  modeling  turbulent  flows  of 
supercritical  fluids. 
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Figure  6-17a.  Comparison  of  Heat  Flux  Predictions  from  the  Three  Turbulence 
Models  for  Upward  Flow  in  the  2-D  channel  with  Tw  =  320°K 


Figure  6-17b.  Comparison  of  Heat  Flux  Predictions  from  the  Three  Turbulence 
Models  for  Upward  Flow  in  the  2-D  Channel  with  q  =  30KW /M2 


7.  COMPUTATIONAL  RESULTS  FOR  FLOW  AND  HEAT  TRANSFER 


This  section  describes  the  flow  and  heat  transfer  results  obtained  for  supercritical 
fluid  flow  of  COz  over  a  range  of  Reynolds  number  going  from  the  laminar  to  the 
turbulent  regime.  The  effect  of  buoyancy  on  laminar  and  turbulent  heat  transfer  is 
investigated  in  detail.  Detailed  comparisons  with  experimental  data  are  also 
presented.  A  large  amount  of  experimental  data  are  available  (for  flow  of 
supercritical  C02)  and  therefore  this  study  is  of  significant  value  in  terms  of  model 
validation. 

7.1  Boundary  Conditions 

The  geometry  (and  computational  grid)  that  was  considered  for  the  two- 
dimensional  channel  flow  cases  is  shown  in  Figure  7-1.  Only  the  lower  half  of  the 
channel  is  considered  in  the  simulations  with  appropriate  symmetry  boundary 
conditions.  The  half  channel  height  (H)  is  10  mm  and  the  channel  length  (L)  is  150 
mm.  The  channel  pressure  was  set  equal  to  7.58  MPa  which  is  above  the  critical 
pressure  (Pc=  7.38  MPa)  for  C02.  The  inlet  temperature  of  C02  was  set  to  288°K 
which  is  lower  than  the  critical  temperature  (Tc  =  304. 1°K)  for  C02.  The  wall  of  the 
channel  is  heated  (Tw  =  320°K)  so  that  the  fluid  temperature  increases  beyond  the 
critical  point  at  some  intermediate  point  along  the  channel.  These  conditions  are 
representative  of  the  conditions  in  fuel  lines  where  the  fuel  (at  supercritical 
pressure  and  subcritical  temperature)  enters  the  heat  exchanger  and  gets  heated 
beyond  the  critical  point. 

It  is  observed  from  Figure  7-1  that  the  computational  grid  is  clustered  near  the 
heated  wall  in  order  to  resolve  the  steep  variations  in  transport  properties.  Grid 
independence  was  checked  by  progressively  refining  the  grid  until  the  solutions 
became  independent  of  the  grid.  Simulations  are  done  for  downward  flows  (in  the 
direction  of  gravity)  and  upward  flows  (against  the  direction  of  gravity)  and  the  heat 
transfer  rates  are  compared  with  the  no-gravity  case.  This  was  done  primarily  to 
assess  the  effort  of  buoyancy  in  influencing  heat  transfer  process. 
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Inlet: 


T  =  288 K  Tc  =  304.1  °K  Channel  Pressure  =  7.58  MPa 

U  =  Laminar_Prof 

Re  =  U  •  2H 


Pc  =  7.38  MPa 


Figure  7-1.  Computational  Geometry  and  Grid  for  the  2-D  Channel  Flow  Problem 


7.2  Results  for  Laminar  Heat  Transfer 

The  inlet  velocity  of  C02  was  taken  to  be  0.05m/ sec.  This  corresponds  to  a  Reynolds 
number  of  10,000  based  on  inlet  conditions.  A  laminar  flow  model  was  used  to 
calculate  the  fluid  flow  and  heat  transfer  in  the  channel.  Upward  and  downward 
flows  and  the  corresponding  no-gravity  cases  were  simulated.  Figure  7-2  shows  the 
transverse  profiles  of  density,  velocity  and  temperature  at  an  axial  location 
corresponding  to  x/2H  =  4.4.  The  dashed  lines  represent  the  no-gravity  case.  The 
heated  wall  causes  the  local  fluid  temperature  at  the  wall  to  increase  beyond  the 
critical  point  although  the  bulk  fluid  is  still  at  subcritical  temperatures.  This  gives 
rise  to  a  density  gradient  near  the  wall.  The  density  decreases  by  almost  a  factor  of 
four  in  the  density  boundary  layer.  For  the  no-gravity  case,  the  density  and 
temperature  boundary  layer  thicknesses  are  of  the  order  of  0.25  in  non-dimensional 
distance.  Interestingly,  the  upward  flow  shows  much  higher  gradients  for  the  same 
quantities.  The  boundary  layer  thicknesses  (for  the  upward  flow)  for  the  density  and 
temperature  are  of  the  order  of  0.05  in  non-dimensional  distance. 

This  feature  can  be  better  understood  by  looking  at  the  velocity  profiles.  The  no 
gravity  cases  exhibits  the  typical  boundary  layer  profile  with  the  velocity  increasing 
monotonically  from  a  zero  value  at  the  wall  to  the  bulk  value  outside  the  boundary 
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layer.  However,  for  the  upward  flow,  the  velocity  profile  shows  a  peak  near  the 
wall.  Also  the  magnitude  of  the  peak  is  quite  large  compared  to  the  bulk  value. 
This  is  due  to  the  fact  that  the  large  decrease  in  density  near  the  wall  gives  rise  to  an 
upward  acceleration  of  the  fluid  under  the  action  of  gravity.  The  local  acceleration 
of  the  low  density  fluid  leads  to  thinner  boundary  layers  for  both  temperature  and 
density.  The  direct  result  of  a  thinner  boundary  layer  for  the  temperature  is  an 
increased  heat  transfer  coefficient.  This  is  reflected  in  Figure  7-4  which  shows  the 
local  heat  flux  at  the  wall  along  the  entire  channel  length.  It  is  observed  that  the 
heat  flux  values  calculated  for  the  upward  flow  are  almost  an  order  of  magnitude 
higher  than  those  computed  for  the  no-gravity  case. 

Figure  7-3  shows  the  transverse  profiles  of  density,  temperature  and  velocity  for  the 
downward  flow  at  the  same  axial  location.  The  results  for  the  no-gravity  case  are 
shown  by  dashed  lines.  Again,  it  is  observed  that  the  effect  of  buoyancy  is  to 
accelerate  the  low  density  fluid  near  the  wall  in  the  upward  direction  (i.e.,  against 
the  direction  of  the  bulk  flow).  The  density  and  temperature  boundary  layer 
thicknesses  for  the  downward  flow  case  are  thinner  compared  to  the  no-gravity  case. 
This  translates  into  higher  heat  transfer  rates  as  shown  in  Figure  7-4.  The  heat 
transfer  rates  at  the  inlet  (for  the  downward  flow)  are  lower  than  those  obtained  for 
the  upward  flow.  The  two  curves  cross  each  other  half-way  along  the  length  of  the 
channel.  This  can  be  explained  by  the  fact  that  the  upward  acceleration  due  to 
buoyancy  counteracts  the  mean  downward  flow  at  the  inlet  leading  to  local  flow 
stagnation  and  decreased  heat  transfer  rates.  However,  at  locations  far  downstream 
of  the  inlet,  the  upward  flow  due  to  buoyancy  overwhelms  the  effect  of  the  bulk 
flow  leading  to  increased  heat  transfer  rates. 
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Figure  7-2.  Predicted  Profiles  of  Density,  Velocity  and  Temperature  in  Channel  for 
the  Upward  Flow  (Re  =  1,000) 


Figure  7-3.  Predicted  Profiles  of  Density,  Velocity  and  Temperature  in  the  Channel 
for  the  Downward  Flow  (Re  =  10,000) 


Figure  7-4.  Local  Heat  Flux  Along  the  Channel  Wall  for  the  Upward  Flow,  the 
Downward  Flow  and  the  No-Gravity  Case  (Re  =  10,000) 

7.3  Discussion  of  Laminar  Heat  Transfer  Results 

It  is  observed  that  the  effect  of  buoyancy  is  to  increase  local  heat  transfer  rates  for 
both  upward  and  downward  flows.  This  is  mainly  due  to  the  local  upward 
acceleration  of  the  low  density  fluid  near  the  heated  wall  under  the  action  of  gravity. 
Since  the  heat  transfer  rates  in  this  flow  regime  are  largely  influenced  by  the  local 
density  gradients,  it  is  a  fair  assumption  that  ideal  gas  models  would  not  be  able  to 
predict  heat  transfer  in  such  flows  mainly  due  to  the  inability  of  ideal  gas  models  in 
predicting  the  large  density  decrease  near  the  heated  wall. 

Figure  7-5  shows  the  comparison  between  predicted  heat  transfer  rates  and 
measured  values  for  laminar  free  convection  of  supercritical  C02  adjacent  to  a 
heated  wall.  The  surrounding  pressure  is  8.11  MPa  and  the  ambient  temperature  is 
310.16°K.  Heat  transfer  measurements  [32]  were  made  while  the  flow  was  still  in  the 
laminar  regime.  The  heat  flux  is  plotted  as  a  function  of  wall  temperature  in  Figure 
7-5.  It  is  observed  that  the  model  predictions  are  in  good  agreement  with  the  data. 


90 


4241/6 


Figure  7-5.  Comparison  of  Measured  and  Predicted  Heat  Flux  Values  for  Free 
Convection  of  Supercritical  C02 

7.4  Results  for  Turbulent  Heat  Transfer 

The  turbulent  flow  simulations  were  done  using  the  low  Reynolds  number  k-e 
model  [20].  The  case  of  laminar  flow  in  the  two-dimensional  channel  was  repeated 
using  the  turbulent  flow  model.  The  heat  transfer  rates  predicted  by  the  laminar 
model  and  the  turbulent  model  are  compared  for  the  upward  flow,  the  downward 
flow  and  the  no-gravity  case  (see  Figure  7-6).  It  is  observed  that  both  models 
essentially  predict  the  same  heat  transfer  values  thus  indicating  that  the  effect  of 
turbulence  was  minimal  for  this  problem. 

The  inlet  velocity  was  increased  by  a  factor  of  ten,  i.e.,  the  effective  Reynolds 
number  was  increased  to  100,000.  Figure  7-7  shows  the  turbulent  heat  transfer  rates 
predicted  along  the  channel  wall  (for  a  wall  temperature  of  320°K)  for  the  upward 
flow,  the  downward  flow  and  the  no-gravity  case.  In  contrast  to  the  laminar  case,  it 
is  observed  that  the  heat  transfer  rates  for  the  upward  flow  are  lower  than  those  for 
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the  no-gravity  case  which  are  in  turn  lower  than  the  heat  transfer  rates  predicted  for 
the  downward  flow.  These  findings  are  contrary  to  the  results  obtained  for  Re  = 
10,000  where  the  predicted  heat  flux  values  for  the  upward  flow  and  the  downward 
flow  are  higher  than  those  for  the  no-gravity  case.  This  rather  curious  behavior  can 
be  explained  based  on  the  following  arguments. 

Figures  7- 8a  and  7-8b  show  the  velocity  profiles  at  different  locations  in  the  channel 
for  the  upward  and  downward  flow,  respectively.  The  dashed  lines  indicate  the 
results  for  the  no-gravity  case.  For  the  upward  flow,  it  is  observed  that  the  upward 
acceleration  of  the  low  density  fluid  near  the  wall  due  to  gravity  gives  rise  to 
velocity  peaks  or  maxima  near  the  wall.  The  net  effect  of  the  velocity  maxima  is  to 
increase  the  velocity  gradient  at  the  wall  and  to  decrease  the  velocity  gradient 
(compared  to  the  no-gravity  case)  a  little  distance  away  from  the  wall.  In  turbulent 
boundary  layers,  a  shear  layer  exists  in  the  vicinity  of  the  wall  that  is  essentially 
responsible  for  the  turbulent  transport  of  mass  momentum  and  energy  from  the 
bulk  flow  to  the  boundary  layer.  But  the  occurrence  of  velocity  maxima  reduces  the 
velocity  gradient  in  the  region  of  the  shear  layer  thus  reducing  the  turbulence 
kinetic  energy  and  hence  the  turbulent  transport.  For  the  downward  flow  (Figure,  7- 
8b)  it  is  observed  that  the  velocity  gradients  in  the  shear  layer  region  are  actually 
increased  (compared  to  the  no-gravity  case)  due  to  the  acceleration  of  the  low 
density  fluid  and  thus  leading  to  increased  turbulent  transport. 
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Figure  7- 8a.  Comparison  of  Mean  Velocity  Profiles  for  the  Upward  Flow  and  the 
No-Gravity  Cases  for  Re  =  100,000 


Figure  7-8b.  Comparison  of  Mean  Velocity  Profiles  for  the  Downward  Flow  and 
No-Gravity  Cases  for  Re  =  100,000 


Figures  7- 9a  and  7- 9b  show  the  turbulent  viscosity  profiles  (calculated  by  the  low 
Reynolds  number  k-e  model)  for  the  upward  and  the  downward  flows  respectively. 
For  the  upward  flow,  it  is  observed  that  the  calculated  turbulent  viscosity  near  the 
wall  (in  the  shear  layer  region)  is  lower  than  those  predicted  for  the  no-gravity  case. 
Correspondingly,  the  turbulent  viscosity  for  the  downward  flow  is  higher  than 
those  obtained  for  the  no-gravity  case.  Figures  7-10a  and  7-1  Ob  show  the 
temperature  profiles  in  the  channel  for  the  upward  flow  and  the  downward  flow 
respectively.  For  the  upward  flow,  it  is  observed  that  the  temperature  gradients  at 
the  wall  are  lower  than  those  predicted  for  the  no-gravity  case.  This  is  a  direct  result 
of  reduced  turbulent  transport  (due  to  velocity  maxima)  and  leads  to  reduced  heat 
transfer  rates  or  increased  heat  transfer  resistance  at  the  wall.  Correspondingly,  the 
temperature  gradients  for  the  downward  flow  are  larger  than  the  gradients  for  the 
no-gravity  case  thus  leading  to  increased  heat  transfer  at  the  wall. 

7.5  Discussion  of  Turbulent  Heat  Transfer  Results 

It  is  observed  that  the  upward  acceleration  of  low  density  fluid  near  the  heated  wall 
(due  to  gravity)  affects  the  local  velocity  profile  thus  influencing  the  turbulent 
transport  near  the  wall.  The  trends  in  heat  transfer  are  also  changed  depending  on 
the  direction  of  the  mean  flow  with  respect  to  gravity.  Similar  phenomena  have 
been  observed  in  experiments  [4]  involving  the  turbulent  flow  of  C02  in  pipes. 

Figures  7-lla  and  7-llb  show  the  results  for  heat  transfer  calculated  in  a  pipe  (1.9  cm 
in  diameter)  with  a  specified  heat  flux  condition.  The  Reynolds  number  for  these 
cases  are  approximately  equal  to  100,000.  Figure  7-lla  shows  the  computed  wall 
temperatures  predicted  for  a  heat  flux  condition  of  3.09  W/cm2.  The  symbols 
represent  experimental  measurements  [4]  with  the  downward  triangles  representing 
downflow  data  and  the  upward  triangles  representing  the  upflow  data.  For  a  heat 
flux  value  of  3.09  W/cm2,  both  the  data  and  the  model  indicate  that  the  wall 
temperatures  for  the  upward  flow  are  slightly  higher  than  those  for  the  downward 
flow.  Again,  this  is  a  direct  consequence  of  increased  heat  transfer  resistance  for  the 
upward  flow.  The  model  results  appear  to  be  in  reasonable  agreement  with  the  data. 
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Figure  7-9 a.  Comparison  of  Turbulent  Viscosity  Profiles  for  the  Upward  Flow  and 
No-Gravity  Cases  for  Re  =  100,000 


Figure  7-9b.  Comparison  of  Turbulent  Viscosity  Profiles  for  the  Downward  Flow 
and  No-Gravity  Cases  for  Re  =  100,000 


f] 


Figure  7-10a.  Comparison  of  Temperature  Profiles  for  the  Upward  Flow  and  No- 
Gravity  Cases  for  Re  =  100,000 


Figure  7-10b.  Comparison  of  Temperature  Profile  for  the  Downward  Flow  and  No- 
Gravity  Cases  for  Re  =  100,000 
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Figure  7- lib  shows  the  results  for  a  higher  heat  flux  value  of  5.19  W/cm2.  The 
predicted  wall  temperatures  for  the  upward  flow  show  a  dramatic  increase 
compared  to  the  downward  flow.  This  is  reflected  by  both  the  model  and  the  data. 
Although  the  trends  are  correctly  predicted  by  the  model,  there  are  some  differences 
in  the  actual  magnitude  of  wall  temperatures  between  the  model  and  the  data.  The 
reason  for  this  discrepancy  is  not  clear.  However,  the  experimental  data  itself  has 
significant  error  bars  (or  standard  error)  for  the  high  heat  flux  conditions  [4] 
reflecting  the  difficulty  of  making  accurate  measurements  in  this  regime.  But  there 
is  no  mistaking  the  trends  portrayed  by  both  the  data  and  the  model  predictions. 

For  the  low  heat  flux  case,  the  heat  addition  is  not  sufficient  to  create  a  significantly 
thick  low  density  layer  near  the  wall  (by  increasing  the  fluid  temperature  beyond  its 
critical  point)  and  hence  only  small  differences  are  observed  between  the  upflow 
and  downflow  results.  However,  in  the  high  heat  flux  case,  the  buoyancy  driven 
flow  is  of  sufficient  magnitude  to  drastically  increase  the  heat  transfer  resistance  for 
the  upflow  case  thus  resulting  in  increased  wall  temperatures  and  local  hot  spots  in 
the  tube. 

7.6  Non-Dimensional  Analysis  of  the  Problem 

One  of  the  main  observations  from  the  above  simulations  is  that  buoyancy  plays  an 
important  role  in  influencing  heat  transfer  (in  supercritical  flows)  even  at  Reynolds 
numbers  as  high  as  100,000.  This  is  of  particular  concern  for  aircraft  heat  exchangers 
which  operate  in  this  range  of  Reynolds  numbers.  In  order  to  better  understand  the 
influence  of  buoyancy,  the  following  non-dimensional  analysis  is  performed. 

The  Grashof  number  is  a  non-dimensional  number  that  is  a  ratio  of  the  effect  of 
buoyancy  forces  over  viscous  forces  in  a  flow.  The  Grashof  number  is  defined  as 
follows: 


^  _  PgApL3 
r"  H2 

where  L  is  a  typical  length  scale  for  the  problem.  For  the  operating  conditions 
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described  for  the  two-dimensional  channel  flow  and  the  pipe  flow,  a  representative 
Grashof  number  can  be  calculated  as  follows  based  on  inlet  conditions: 

p  =  860  kg/  m3(at  288°K) 
g  =  9.81  m/sec  2 

H  ~  1  x  10-4  N  -  sec  /m2 
Ap  =  p(288°K)  -  pw(320°K) 

*  860  kg/m3  -  200kg/ m3  =  660kg/m3 
L  =  0.15  m 


The  calculated  Grashof  number  turns  out  to  be: 

Gr=  1012 

The  non-dimensional  ratio  that  describes  the  effectiveness  of  natural  convection 
over  forced  convection  is: 


Gr 

(Re)2 


For  a  Reynolds  number  value  of  100,000,  the  ratio  is  still  of  the  order  of  100 
indicating  that  buoyancy  will  play  an  important  role.  This  has  been  confirmed  by 
detailed  simulations  of  the  channel  flow  and  pipe  flow  cases. 

However,  for  a  Reynolds  number  of  1,000,000,  the  above  ratio  will  be  of  the  order  of 
unity  and  the  effect  of  buoyancy  should  be  considerably  smaller.  In  order  to  verify 
this,  additional  simulations  were  performed  (for  wall  temperatures  of  320°K  and 
600  K)  for  the  two-dimensional  channel  with  increased  inlet  velocity  corresponding 
to  a  Reynolds  number  of  1,000,000.  Figure  7-12  shows  the  calculated  heat  transfer 
rates  for  the  two  wall  temperature  cases.  It  is  observed  that  the  heat  transfer  rates 
predicted  for  the  upward  flow  are  exactly  the  same  as  those  predicted  for  the  no¬ 
gravity  case,  i.e.,  gravitational  effects  are  negligible. 
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For  most  cases,  however,  the  evaluation  of  a  representative  Grashof  number 
becomes  extremely  difficult  since  the  density,  the  density  gradient  and  the  viscosity 
become  complex  functions  of  the  reduced  pressure  and  temperature,  the  g  value  and 
flow  characteristics.  Especially  for  problems  with  specified  heat  flux  conditions 
(which  is  representative  of  many  practical  problems  as  opposed  to  the  constant  wall 
temperature  condition),  a  non-dimensional  analysis  of  the  kind  shown  above  is  no 
trivial  task.  The  local  density  gradient  influences  the  flow  field  (through  buoyancy 
effects)  which  in  turn  influences  heat  transfer  (by  affecting  boundary  layers  and  heat 
transfer  coefficients)  which  in  turn  influences  the  local  temperature  distribution 
which  in  turn  affects  the  local  density  distribution.  This  is  illustrated  schematically 
in  Figure  7-13.  The  link  between  temperature  and  density  (as  illustrated  in  Figure  7- 
13)  is  crucial.  For  incompressible  flows  with  low  heat  transfer  rates,  the  link  is  very 
weak.  For  low  speed  flows  of  gases  (in  the  range  of  pressures  where  the  ideal  gas 
assumption  is  valid)  the  density  is  proportional  to  the  reciprocal  of  the  temperature. 
For  supercritical  flows,  the  density  is  a  highly  non-linear  function  of  temperature  in 
the  vicinity  of  the  critical  point. 


Figure  7-12.  Predicted  Heat  Flux  Values  Along  the  Channel  Wall  for  Re  =  1,000,000 
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Figure  7-13.  Illustration  of  the  Coupling  between  Various  Physical  Phenomena 

This  highly  non-linear  coupling  between  temperature  and  density  (in  the 
supercritical  regime)  makes  the  a  priori  determination  of  buoyancy  effects  extremely 
difficult.  If  the  heat  flux  is  sufficient  to  increase  the  fluid  temperature  beyond  the 
critical  value,  the  corresponding  Grashof  numbers  will  be  large.  On  the  other  hand, 
if  the  heat  flux  is  not  strong  enough  to  increase  the  fluid  temperature  beyond  its 
critical  value,  local  Grashof  numbers  will  be  relatively  small.  Therefore  depending 
on  the  heat  flux  value,  Grashof  numbers  can  differ  by  orders  of  magnitude.  For 
example,  if  the  wall  temperature  were  to  be  reduced  from  320°K  to  300°K,  the 
calculated  Grashof  number  is  of  the  order  of  1011,  i.e.,  a  decrease  of  20°K  (~  6%)  in 
temperature  has  reduced  the  Gr  value  by  one  order  of  magnitude.  These  difficulties 
emphasize  the  need  for  advanced  models  that  can  simulate  local  effects  in 
supercritical  flows. 
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8.  ASSESSMENT  OF  THE  THERMAL  STABILITY/DEPOSITION  MODEL 


This  section  presents  the  results  from  the  application  of  the  thermal  stability  models 
discussed  in  Section  2.  The  model  predictions  are  compared  with  experimental  data. 
These  comparisons  are  used  to  further  refine  and  recalibrate  the  models  in  order  to 
improve  their  generality. 

8.1  Tet-A  Experiments 

The  Jet- A  experiments  are  described  in  detail  in  the  article  by  Katta  et  al.  [23].  The 
near  isothermal  cases  were  considered  for  the  simulation  because  it  eliminates  any 
uncertainties  arising  from  inaccuracies  in  the  calculation  of  the  temperature.  The 
operating  conditions  of  the  experiments  are  shown  below: 

Fuel  Type:  Jet-A 
Flow  Rate:  0.5  cc/ min 
Inlet  Temperature:  288K 
Tube  Wall  Temperature:  458K 

Inlet  Reynolds  Number:  Re  =  )  =  2 

Two  tube  sizes  were  considered  for  this  study,  the  first  one  (the  standard  size)  has  an 
inner  diameter  of  2.159  mm  and  the  other  (the  large  size)  has  an  inner  diameter  of 
4.32  mm.  The  9-step  Katta  model  (described  in  Section  2)  was  applied  to  simulate 
these  experiments.  Figures  8-la  and  8-lb  show  the  results  for  the  deposition  rate 
and  the  species  distributions  in  the  tubes.  The  dashed  lines  indicate  the  results 
predicted  by  the  9-step  Katta  model  incorporated  in  CFD-ACE.  Figure  8-2  shows  the 
results  from  the  article  by  Katta  et  al.  [23].  It  is  observed  that  the  CFD-ACE 
predictions  are  in  good  agreement  with  the  results  reported  by  Katta  et  al. 
However,  it  is  observed  that  the  model  is  not  able  to  predict  correctly  the  location  of 
the  deposition  peaks  compared  with  experimental  data.  This  is  a  characteristic 
feature  of  the  9-step  model  that  is  prevalent  in  all  of  the  cases  analyzed  during  this 
study. 

Although  the  results  are  in  good  agreement  with  the  values  reported  in  the 
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publication  by  Katta  et  al.  [23],  the  current  model  shows  small  discrepancies  from 
Katta's  results.  These  deficiencies  were  attributed  to  the  sensitivity  of  the  Katta 

model  to  the  value  of  xw  used  in  the  deposition  model  (Equation  (48)  in  Section  2). 
In  the  original  calculations  reported  in  Katta  et  al.,  a  standard  k-g  model  was  used  for 
the  simulations  and  xw  was  based  on  the  velocity  of  the  first  cell  from  the  wall. 
However  for  large  values  of  xw,  (i.e.  for  turbulent  flows)  the  model  accuracy 
deteriorates.  Katta  et  al.  corrected  for  this  effect  by  placing  a  limiter  on  xw.  However, 
this  procedure  limits  the  generality  of  the  model. 

This  is  demonstrated  by  applying  the  model  to  a  fuel  flow  experiment  described  in 
Section  3.  The  fuel  flows  through  a  heated  tube  with  a  wall  temperature 
distribution  as  shown  in  Figure  8-3b.  The  model  predictions  are  compared  with 
measured  data  (see  Figure  8-3a)  on  the  deposition  rate.  It  is  observed  that  the  model 
(see  Figure  8-3a)  underpredicts  the  total  deposition  rate  as  well  as  the  distribution. 
In  order  to  improve  its  generality,  the  surface  reactions  in  the  Katta  model  were 
recalibrated  (see  Section  2).  The  results  of  the  recalibrated  model  are  shown  by  the 
dark  curves  in  Figures  8-1  and  8-3a.  It  is  observed  that  the  recalibrated  model  yields 
about  the  same  results  for  the  lowspeed  isothermal  tube  cases  (Figures  8-la  and  8-lb) 
as  well  as  better  agreement  in  deposition  rate  for  the  high  speed  fuel  case  shown  in 
Figure  8-3a.  However,  the  location  of  the  deposition  peaks  are  not  predicted 
correctly.  This  is  to  be  expected  since  the  original  9-step  model  shows  the  same 
behavior.  All  of  the  Jet-A  deposition  studies  reported  in  the  subsequent  sections  of 
the  report  were  done  using  the  recalibrated  Katta  model. 

8.2  TP-5  Experiments 

The  details  of  the  JP-5  experiment  conducted  by  Marteney  and  Spadaccini  [33]  are 
described  below: 

Fuel  Type:  JP-5 
Flow  Rate:  2.13  m/s 
Inlet  Temperature:  300°K 
Wall  Heat  Flux:  395  kw/m2 
Inlet  Reynolds  Number:  ~3,000 
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The  flow  of  JP-5  in  a  channel  with  a  specified  heat  flux  at  the  walls  was  considered. 

The  low  Reynolds  number  k-e  model  was  used  to  simulate  the  turbulence  in  the 
flow  field.  The  9-step  model  (for  Jet-A)  as  well  as  the  4-step  Krazinski  model  for  JP-5 
(see  Section  2)  were  applied  to  simulate  this  case.  The  results  are  shown  in  Figures 
8-4a  and  8-4b.  Figure  8-4a  shows  the  wall  and  bulk  temperature  distributions  in  the 
tube.  It  is  observed  that  the  temperature  distribution  is  predicted  reasonably  well  by 

the  low  Re  k-e  model.  The  standard  k-e  model  is  unable  to  predict  the  initial  bump 
in  the  wall  temperature  distribution  due  to  the  inherent  assumption  of  a  fully 
developed  turbulent  flow  profile  at  the  inlet. 

The  results  for  deposition  rate  indicate  that  the  4-step  model  for  JP-5  shows  good 
agreement  with  the  data.  The  9-step  model  (for  Jet-A)  underpredicts  the  deposition 
rate.  The  surface  reactions  in  the  9-step  model  were  further  recalibrated  (see  Section 
2)  for  JP-5  in  order  to  improve  the  agreement.  The  recalibrated  9-step  model  for  JP-5 
shows  better  agreement  with  data.  It  is  expected  that  further  recalibration  of  gas 
phase  reactions  for  JP-5  will  result  in  better  agreement  with  data. 

8.3  Conclusions 


It  was  determined  that  the  original  9-step  Katta  model  had  to  be  recalibrated  to 
improve  its  generality  especially  for  higher  fuel  flow  rates.  The  recalibrated  model 
showed  good  agreement  with  the  Katta  model  for  low  speed  flows  as  well  as  for 
high  flow  rate  fuel  experiments.  The  9-step  model  needed  further  recalibration  for 
JP-5  in  order  to  match  the  deposition  rate  data. 
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Deposition  (pg/cm2/h) 


Axial  Distance  (cm) 

Figure  8-2.  Predicted  and  Measured  Deposition  in  Near-Isothermal  (458K)  Heated- 
Tube  Experiment  Using  Two  Different  Bore  Sizes  [23] 
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Comparison  of  Deposition  Rates  for  Jet-A  Case  7/94-6 
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Figure  8-3a.  Deposition  Rate  Distribution  Along  the  Tube 

Comparison  of  Wall  Temperatures  for  Jet-A  Case  7/94-6 


Figure  8-3b.  Wall  Temperature  Distributions  Along  the  Tube 


Temperature  Distribution  of  JP5  Case 


Axial  Distance  (m) 

Figure  8-4a.  Wall  and  Bulk  Temperature  Distributions  Along  the  Tube 


Deposition  Mechanisms  for  JP-5  Case 


1000/Tw,tl  (1/K> 

Figure  8-4b.  Deposition  Rate  as  a  Function  of  Wall  Temperature  for  Different 
Chemistry  Models 
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9. 


BENCHMARK  VALIDATION  STUDY  FOR  THE  DEPOSITION  MODEL 


This  section  presents  the  comparisons  between  model  predictions  and  data  obtained 
at  Wright  Laboratory  during  the  summer  months  (June  to  August)  of  1994.  Heat 
transfer  and  deposition  experiments  were  conducted  with  air  saturated  Jet-A  2926  as 
the  test  fuel.  The  details  of  the  experimental  conditions  and  apparatus  are  described 
in  Section  3.  Three  test  cases  (covering  two  operating  pressures  and  two  flow  rates) 
were  chosen  (from  all  the  sets  of  experiments)  as  benchmark  cases  for  the  validation 
study.  These  cases  are  as  follows: 

a.  Case  7/94-6  (Pressure  =  650  psia,  flow  rate  =  102  cc/ min) 

b.  Case  7/94-10  (Pressure  =  400  psia,  flow  rate  =  105  cc/ min) 

c.  Case  7/94-2  (Pressure  =  650  psia,  flow  rate  =  190  cc/min) 

These  were  short  duration  experiments  of  5  hours.  Numerical  simulations  of  the 
above  experiments  were  performed  and  the  model  predictions  of  outer  wall 
temperature  and  fuel  deposition  rate  were  compared  with  measured  data. 

9.1  Description  of  the  Computational  Model 

The  schematic  of  the  computational  domain  is  shown  in  Figure  9-1.  The  top  surface 
is  the  hot  furnace  wall.  The  space  between  the  furnace  wall  and  the  steel  tube 
consists  of  stationary  air.  The  heat  is  transferred  from  the  furnace  to  the  steel  tube 
by  means  of  radiation  as  well  as  conduction  through  the  air.  The  steel  tube  is  of 
finite  thickness  and  offers  conduction  resistance  to  the  flow  of  heat.  The  fuel  flow 
in  the  steel  tube  is  turbulent. 
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Figure  9-1.  Schematic  of  Furnace  and  Tube  Geometry 

The  furnace  temperature  is  set  to  either  1200°C  or  1000°C  depending  upon  the 
experiment  to  be  modeled.  The  thermal  conductivity  and  specific  heat  of  the 
stationary  air  are  taken  to  be  0.06  W/m-k  and  1000  J/kg-K,  respectively.  The  heat 
transfer  inside  the  furnace  is  simulated  using  radiation  and  thermal  conduction, 
with  absorptivity  of  air  as  0.1  and  emissivities  of  furnace  inner  surfaces  as  0.8. 
Radiative  and  conjugate  heat  transfer  through  the  steel  tube  are  modeled.  The 
properties  of  the  tube  are  taken  to  be  that  of  AISI  316  steel  at  650°K,  i.e,  the  thermal 
conductivity  and  specific  heat  are  taken  as  20  W/m-K  and  550  J/Kg-K,  respectively. 
The  laminar  transport  properties  of  the  Jet-A  fuel  are  assumed  to  be  the  same  as 
those  of  n-C12H23.  The  low  Reynolds  number  k-e  model  is  used  to  compute  the 

turbulent  transport.  The  turbulent  conductivity  and  diffusivity  of  the  flow  are 
calculated  from  the  turbulent  viscosity  by  assuming  a  constant  turbulent  Prandtl 
number  of  0.93  and  a  turbulent  Schmidt  number  of  0.9.  The  deposition  rate  is 
calculated  from  the  9-step  chemistry  mechanism  proposed  by  Katta  et  al  (described 
in  Section  2)  with  appropriate  recalibration  of  the  surface  reaction  steps. 
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9.2  Model  Predictions  and  Comparisons  with  Data 


A  grid  independence  study  was  performed  to  assess  the  sensitivity  to  grid 
refinement.  Two  grid  sizes  (75x64  and  120x90)  were  considered.  The  finer  grid  had  a 
greater  clustering  of  points  in  the  vicinity  of  the  temperature  peak.  The  predicted 
wall  temperatures  and  the  deposition  rates  of  the  two  grids  are  almost  the  same 
indicating  that  grid  independence  has  been  achieved  (Figures  9-2a  and  9-2b).  The 
calculations  described  hereafter  have  been  done  on  the  75x64  grid. 

Case  7/94-6:  The  calculated  wall  temperatures  and  deposition  rates  for  this  case  are 
shown  in  Figures  9-3a  and  9-3b,  respectively.  The  trend  and  order  of  magnitude  of 
the  predictions  compare  well  with  experimental  data.  However,  the  predicted 
temperature  and  deposition  peak  are  sharper  and  occur  further  upstream  (close  to 
the  furnace  collar)  than  the  measured  data.  There  may  be  a  number  of  reasons  for 
this  discrepancy  such  as  (i)  the  difficulty  in  modeling  the  boundary  conditions  for 
energy  accurately  near  the  ends  of  the  furnace,  (ii)  the  current  axisymmetric 
assumption  precludes  the  consideration  of  3-dimensional  turbulence  effects  in  the 
tube,  and  (iii)  the  validity  of  the  deposition  model  at  these  temperatures  may  be  in 
question  since  the  calibration  was  done  for  much  lower  temperatures.  However, 
there  is  good  qualitative  and  to  some  extent,  good  quantitative  agreement  between 
the  model  predictions  and  data. 

Case  7/94-10:  This  case  is  similar  to  the  earlier  one  except  that  the  operating 
pressure  is  reduced  to  400  psia.  The  results  are  shown  in  Figures  9-4a  and  9-4b.  The 
wall  temperature  predictions  are  similar  to  results  for  the  earlier  case.  The  peak  in 
deposition  rate  is  still  predicted  further  upstream  than  the  measured  data. 
However,  the  magnitudes  of  the  temperatures  and  deposition  rates  agree  well  with 
the  data. 

Case  7/94-2:  For  this  case,  the  flow  rate  is  increased  to  190  cc/min.  The  results  are 
shown  in  Figures  9-5a  and  9-5b.  Again,  the  results  are  in  qualitative  agreement  with 
data.  The  temperature  and  deposition  rate  peaks  are  predicted  to  be  further 
upstream  as  compared  with  data. 
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Figure  9-3a.  Comparison  of  Predicted  Wall  Temperatures  for  Case  7/94-6  with  Data 
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Figure  9-3b.  Predicted  Deposition  Rates  for  Case  7/94-6  with  Data 
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Figure  9-4a.  Predicted  Wall  Temperatures  for  Case  7/94-10 
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Figure  9-4b.  Predicted  Deposition  Rates  of  Case  7/94-10 
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Figure  9-5a.  Predicted  Wall  Temperatures  for  Case  7/94-2 
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Figure  9-5b.  Predicted  Deposition  Rates  for  Case  7/94-2 


9.3  Discussion  of  Results 


It  is  important  to  understand  the  precise  reason  for  the  discrepancies  between 
predictions  and  the  measured  data.  One  significant  problem  has  been  the 
determination  of  the  accuracy  of  the  measurements.  In  some  instances  the 
measured  values  of  wall  temperature  exceeded  the  furnace  temperature  which  is 
physically  unrealistic.  However,  the  models  also  need  further  improvement  in  the 
following  areas:  (i)  three-dimensional  effects  due  to  buoyancy  induced  turbulence 
need  to  be  considered,  (ii)  the  chemistry  model  needs  to  be  recalibrated  for  the  high 
temperature  regime,  and  (iii)  the  original  Katta  model  did  not  correctly  predict  the 
location  of  the  deposition  peaks  (see  Figure  8-2)  which  indicates  that  further 
improvement  in  the  chemistry  model  is  necessary. 

Analyses  of  the  experimental  data  show  multiple  peaks  in  deposition  rate  along  the 
fuel  tube  during  the  experiments.  This  phenomenon  is  also  predicted  by  the  models 
used  in  the  study.  The  following  discussion  presents  an  explanation  for  the 
multiple  peaks  by  examining  the  local  values  of  temperatures,  oxygen  concentration 
and  turbulent  kinetic  energy. 

The  operating  conditions  are  chosen  to  correspond  to  the  experimental  case  7/94-6 
except  that  a  uniform  heat  flux  is  assumed  at  the  tube  wall.  The  heat  flux  is  treated 
as  a  parameter  in  this  study.  Three  different  simulations  corresponding  to  heat  flux 
values  of  200,  300,  and  400  KW/m2  are  performed.  A  laminar  velocity  profile  is 
assumed  at  the  inlet. 

Figure  9-6a  shows  the  wall  and  bulk  fluid  temperatures  along  the  length  of  the  tube 
for  the  three  cases.  It  is  observed  that  all  cases  show  a  sharp  peak  in  the  wall 
temperature  a  short  distance  downstream  of  the  inlet.  The  reason  for  this  peak  is 
that  the  flow  is  essentially  laminar  in  the  inlet  region  and  the  smaller  diffusion 
transport  causes  the  wall  temperatures  to  increase  sharply  (see  also  Figures  9-7a,  9- 
8a,  and  9-9a).  However,  as  the  flow  becomes  turbulent  (see  turbulent  kinetic  energy 
distributions  in  Figures  9-7c,  9-8c,  and  9-9c)  the  diffusion  transport  increases  and  the 
wall  temperatures  show  an  abrupt  decrease  from  the  peak  value  and  a  gradual 
increase  along  the  tube.  Figure  9-6b  shows  the  deposition  rate  (on  the  wall)  along 
the  tube  for  the  three  cases.  It  is  observed  that  there  is  a  peak  close  to  the  inlet 
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followed  by  another  peak  at  a  downstream  location.  Comparing  Figures  9-6a  and  9- 
6b  it  is  obvious  that  the  first  peak  is  caused  by  the  sudden  rise  in  the  temperature 
near  the  inlet.  The  local  hot  spot  will  result  in  the  reaction  of  the  dissolved  02  (see 
Figures  9-7b,  9-8b  and  9-9b)  to  form  the  deposition  precursor  very  near  the  wall  to 
give  rise  to  the  first  deposition  peak. 

The  formation  of  the  second  peak  however,  is  more  dominated  by  what  happens  in 
the  bulk  flow,  i.e.,  it  is  not  a  local  process.  As  the  flow  becomes  turbulent  and  the 
bulk  fluid  gets  heated,  the  reaction  in  the  bulk  fluid  gives  rise  to  precursors.  These 
precursors  are  transported  to  the  wall  by  turbulent  diffusion  giving  rise  to  the 
second  peak.  Figures  9-7b,  9-8b  and  9-9b  illustrate  the  occurrence  of  the  deposition 
peaks  by  looking  at  gradients  in  oxygen  concentration.  It  is  obvious  from  the  figures 
that  the  first  peak  is  a  local  phenomena  while  the  second  peak  corresponds  to 
oxygen  concentration  gradients  in  the  bulk  fluid. 

Figures  9-10a  and  9-10b  show  the  computed  results  for  a  heat  flux  of  300  kW/m2  for 
three  different  surface  reaction  models.  The  first  model  considers  only  reaction  7  on 
the  surface,  the  second  model  considers  only  reactions  8  and  9  on  the  surface  and  the 
third  model  considers  reactions  7,  8  and  9.  It  is  observed  from  Figure  9-10a  that 
reaction  7  by  itself  is  unable  to  predict  the  deposition  peaks.  The  model  with 
reactions  8  and  9  only  is  very  close  to  the  full  model  indicating  that  these  reactions 
are  the  dominant  mechanisms  for  deposition.  Figure  9-10b  shows  the  deposition 
rates  on  a  log  scale.  It  is  observed  that  at  deposition  rates  smaller  than  1.0|ig/cm2-hr, 
reaction  7  is  dominant.  For  higher  deposition  rates,  reactions  8  and  9  (reaction  8  in 
particular)  are  the  major  contributors  to  the  bulk  production  of  the  precursor. 
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Figure  9-7.  Results  for  qw  =  200  kW  / m2 
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(c)  turbulent  kinetic  energy  distribution 


Figure  9-8.  Results  for  qw  =  300  kW/m2 
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(b)  oxygen  mass  fraction  distribution 


Figure  9-9.  Results  for  qw  =  400  kW/m2 
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Deposition  Rale  (p9/cm2  .^) 


(a)  temperature  distribution 


9-10.  Computational  Results  with  3  Different  Surface  Reaction  Models 
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10.  BENCHMARK  VALIDATION  STUDY  FOR  HEAT  TRANSFER  MODEL 


This  section  presents  the  results  from  a  series  of  simulations  done  to  model  the 
experiments  conducted  at  the  University  of  Iowa  on  heat  transfer  in  SF6  flow. 

10.1  SFg  Experiments 

In  the  experiment,  pressurized  SF6  enters  a  thin  stainless  steel  tube,  which  is  heated 
by  radiation  from  a  furnace.  Test  conditions  included  two  flow  rates  and  reduced 
pressures  ranging  from  0.75  to  1.5.  The  details  of  the  experiment  are  described  in 
Section  4.  Test  cases  corresponding  to  flow  rates  of  150  and  450  ml/min,  at  pressures 
of  400,  531  and  700  psia  are  simulated.  The  predictions  of  tube  wall  and  bulk  fluid 
temperature  at  the  flow  exit  are  compared  with  the  measurements. 

The  set-up  of  the  numerical  problem  is  illustrated  in  Figure  10-1.  The  chamber  is 
assumed  to  be  filled  with  stationary  air  and  heat  transfer  is  accomplished  by 
radiation  (zero  absorptivity)  and  thermal  conduction  (conductivity  of  0.05  W/mK 
and  specific  heat  is  1000  J/kgK).  Conjugate  heat  transfer  is  modeled  at  air/tube,  and 
tube/SF6  interfaces.  The  thermal  properties  of  the  tube  are  assumed  as  conductivity 
=  20  W/mK  and  specific  heat  =  550  J/kgK.  The  turbulence  is  modeled  by  the  low 
Reynolds  number  k-e  model  and  a  turbulent  Prandtl  number  of  0.9  is  assumed. 

10.2  Computational  Results 

The  computational  results  and  comparison  with  data  are  described  below  for  the 
two  flow  rates. 

Flow  Rate  150ml/min:  This  flow  rate  corresponds  to  an  inlet  Reynolds  number  of 
about  2,000.  The  results  for  temperature  at  flow  exit  vs.  furnace  temperature  are 
presented  in  Figure  10-2.  From  the  enthalpy  curve  for  SF6  (see  Figure  5-4d),  it  is 
observed  that  a  large  amount  of  energy  input  is  needed  to  get  the  temperature  over 
the  pseudo-critical  point  for  pressures  of  400  and  531  psia.  This  is  reflected  in  the 
bulk  temperature  predictions  (Figures  10-2  and  10-3)  which  are  flat  at  30  and  45°C 
(302  and  318°K),  the  pseudocritical  temperatures  for  400  and  531  psia,  respectively. 
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However,  the  measured  bulk  temperatures  show  an  increasing  trend.  For  a 
pressure  of  700  psia  (see  Figure  10-4),  the  enthalpy  curve  shows  a  continuous 
increase,  and  the  predicted  bulk  temperature  is  increasing  steadily,  but  still  lower 
than  the  measured  value. 

For  wall  temperatures,  the  predictions  show  similar  trends  as  the  data  (for  400  and 
531  psia),  except  that  the  predicted  wall  temperatures  are  higher  than  the  data.  At  700 
psia,  however,  the  predicted  values  are  lower  than  the  data.  This  may  be  accounted 
for  by  the  trend  of  increasing  laminar  thermal  conductivity  with  increase  in 
pressure. 

Flow  Rate  450  ml/min;  This  flow  rate  corresponds  to  an  inlet  Reynolds  number  of 
about  10,000  and  the  turbulent  transport  effect  is  more  dominant  than  in  the 
previous  cases.  Figures  10-5, 10-6  and  10-7  show  the  results  for  pressures  of  400,  531 
and  700  psia,  respectively.  These  results  also  show  lower  bulk  temperature 
predictions  compared  to  the  data  but  the  increasing  trend  is  captured.  For  the  wall 
temperature,  the  predicted  values  are  lower  than  the  measured  data,  for  furnace 
temperatures  higher  than  600°C.  It  needs  to  be  pointed  out  that  in  the  experiment, 
the  thermo-couples  attached  to  the  steel  tube  acted  as  extra  conductive  wires  for  heat 
transfer  to  the  tube.  Also,  it  was  found  the  temperature  readings  from  thermo¬ 
couples  were  higher  (up  to  50%)  than  the  actual  value,  especially  at  high 
temperatures.  Considering  the  uncontrolled  experimental  conditions  and 
measurement  error,  the  predicted  temperature  curves  are  quantitatively  comparable 
to  the  experimental  results. 
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Figure  10-2.  Results  for  p  =  400  psia  and  Q  =  150  ml/min 
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Figure  10-4.  Results  for  p  =  700  psia  and  Q  =  150  ml/min 


11.  COMPUTATIONAL  RESULTS  FROM  THE  LES  STUDY 


This  section  presents  the  computational  results  obtained  from  the  Large  Eddy 
Simulation  (LES)  study  of  the  flow  of  supercritical  C02  in  a  tube.  The  main 
objective  of  this  study  was  to  resolve  inherent  unsteadiness  (if  any)  in  flows  of 
supercritical  fluids. 

11.1  Problem  Definition 

The  problem  of  C02  flow  in  a  channel  is  considered  (see  schematic  in  Figure  11-la). 
This  problem  has  been  described  in  detail  in  Section  7.  The  pressure  in  the  channel 
is  above  the  critical  pressure  of  C02.  The  inlet  temperature  is  below  the  critical 
temperature.  However,  heating  from  the  constant  temperature  wall  causes  the  fluid 
temperature  to  increase  beyond  the  critical  value.  This  flow  was  studied  in  detail 
(see  Section  7)  for  a  range  of  Reynolds  numbers  and  for  upward  as  well  as 
downward  flow  configurations  (with  respect  to  gravity). 

The  computational  grid  for  the  LES  study  is  shown  in  Figure  11-lb.  The  grid  is 
made  considerably  finer  (compared  to  the  one  used  in  the  earlier  study  described  in 
Section  7)  in  order  to  resolve  any  inherent  unsteadiness  (due  to  small  scale 
structures)  in  these  flows.  The  transient  governing  equations  (i.e.,  with  the  time 
dependent  terms)  are  solved  for  this  problem  with  a  time  step  of  AT  =  0.01  secs.  In 
addition,  the  temperature  is  monitored  at  several  locations  in  the  channel  to 
determine  the  nature  of  time  dependence.  Figure  11-lc  shows  the  monitoring 
locations.  The  locations  1,  2,  and  3  are  placed  very  close  to  the  wall  to  monitor  near 
wall  effects. 

11.2  Results  from  the  LES  Study 

A  series  of  simulations  were  performed  for  three  Reynolds  numbers  of  10,000, 
100,000  and  1,000,000,  respectively.  For  each  Reynolds  number,  three  cases 
corresponding  to  upward  flow  (flow  against  gravity),  downward  flow  (flow  in  the 
direction  of  gravity)  and  horizontal  flow  (flow  perpendicular  to  gravity)  were 
analyzed. 
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Figure  11-2  shows  the  computed  heat  flux  along  the  wall  (for  Re  =  10,000)  for  the 
three  flow  configurations.  The  upward  and  downward  flow  cases  converged  to 
steady  state  solutions  and  the  predicted  values  of  heat  flux  compared  exactly  with 
the  values  obtained  earlier  (see  Section  7)  with  the  steady  state  (coarse  grid)  model. 
However,  the  case  of  horizontal  flow  shows  a  significant  amount  of  flow 
unsteadiness  and  the  computed  heat  flux  (in  Figure  11-2)  shows  large  spatial 
oscillations.  Figure  11-3  shows  the  instantaneous  temperature  distribution  in  the 
channel.  It  is  observed  that  mushroom  shaped  structures  emanating  from  the  wall 
give  rise  to  flow  unsteadiness.  This  has  been  observed  in  experiments  where  the 
large  density  gradient  near  the  wall  (due  to  the  fluid  becoming  supercritical) 
interacts  with  gravity  to  produce  buoyancy  driven  structures.  Figure  11-4  shows  the 
density  distribution  in  the  channel.  It  is  observed  that  the  density  near  the  wall  is 
lower  (by  almost  a  factor  of  four)  than  the  bulk  density  due  to  the  fluid  becoming 
supercritical. 

Figures  ll-5a-c  show  the  calculated  temperatures  as  functions  of  time  at  the  different 
monitoring  locations.  It  is  observed  that  the  buoyancy  driven  flow  causes  small 
scale  eddies  to  be  shed  off  the  wall.  The  monitoring  locations  near  the  wall  (Figure 
ll-5a)  show  a  significant  amount  of  intermittency.  However,  the  temperatures  at 
locations  4,  5  and  6  (0.1  cm  above  the  wall)  shown  in  Figure  ll-5b  indicate  that  the 
frequency  of  intermittency  at  these  locations  is  approximately  half  that  of  locations 
1,  2  and  3.  This  indicates  that  the  small  scale  eddies  near  the  wall  undergo  pairing 
(or  merging)  a  certain  distance  away  from  the  wall  to  give  rise  to  the  lower 
frequencies.  Figure  11 -5c  shows  the  temperature  locations  7,  8  and  9  (0.2  cm  from 
the  wall).  The  frequency  of  intermittency  is  approximately  the  same  as  what  was 
observed  at  locations  4,  5,  and  6. 

The  cases  corresponding  to  Re  =  100,000  and  Re  =  1,000,000  did  not  show  any 
transient  behavior  irrespective  of  the  flow  orientation,  i.e.,  they  converged  to  steady 
state  solutions  that  matched  exactly  with  the  results  obtained  from  the  steady  state 
model. 
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Heat  Flux  (w/m2) 


Figure  11-2.  Predictions  of  Heat  Flux  Along  the  Heated  Wall  of  the  Channel 
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Figure  ll-3b.  Magnified  View  of  Temperature  Distribution  Near  the  Heated  Wall 


134 


4241/6 


DENS  CONTOURS 
FMIN  2.148E+02 


16  7.101  E+02 

18  7.761  E+02 

20  8.422E+02 


Figure  ll-4a.  Density  Distribution  in  the  Channel 
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Figure  ll-5c.  Transient  Temperature  Behavior  at  Locations  (7,  8  and  9)  0.2  cm  Above 
the  Heated  Wall 

11.3  Discussion 

The  results  show  that  the  only  unsteadiness  that  was  observed  in  the  simulations 
was  that  produced  by  buoyancy  effects.  Also,  the  unsteadiness  was  observed  only  at 
low  Reynolds  numbers  (of  the  order  of  10,000)  and  for  horizontal  flow 
configurations.  For  higher  Reynolds  numbers  and  for  upward/downward  flow 
configurations,  the  transient  model  converged  to  steady  state  solutions  that  matched 
the  results  produced  by  the  steady  state  model. 

Therefore,  the  transient  simulations  do  not  show  any  inherent  unsteadiness  in 
supercritical  flows  other  than  buoyancy  effects.  Also,  it  is  clear  that  it  is  not 
necessary  to  invoke  the  transient  model  if  the  steady  state  model  yields  a  converged 
solution.  Only  for  cases  where  the  steady  state  model  fails  to  yield  a  converged 
solution  (i.e.,  for  the  case  of  Re  =  10,000  and  horizontal  flow)  will  the  transient 
model  be  useful  in  terms  of  predicting  intermittency  and  the  local  frequency  of 
unsteady  behavior. 
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12.  CONCLUSIONS 


This  section  presents  the  conclusions  from  the  Phase  II  study  and  plans  for  future 
extensions  and  applications  of  the  models  developed  during  this  study. 

12.1  Summary  of  Phase  II  Work 

The  completion  of  the  Phase  II  study  has  resulted  in  the  following 
accomplishments: 

Incorporation  of  a  Real  Fluid  Transport  Model:  Models  for  computation  of 
thermophysical  properties  of  real  fluids  were  incorporated  into  a  general  purpose 
CFD  code.  These  models  facilitate  the  computation  of  local  thermodynamic  and 
transport  properties  as  a  function  of  pressure,  temperature  and  composition.  Steep 
variations  in  these  properties  at  the  critical  (or  pseudo-critical)  point  can  be  resolved 
reasonably  well  by  this  model.  A  detailed  database  has  been  incorporated  into  the 
code  to  handle  a  range  of  fluids  including  a  number  of  hydrocarbon  fuels.  The  CFD 
code,  CFD-ACE,  is  a  general  purpose  transport  code  that  solves  the  fluid  flow,  heat 
transfer  and  chemistry  in  complex  engineering  systems.  The  turbulence  models  in 
CFD-ACE  were  modified  appropriately  to  account  for  the  flow  of  supercritical  fluids. 

Computation  of  Heat  Transfer  in  Supercritical  Fluids:  The  models  were  applied  to 
analyze  heat  transfer  in  supercritical  fluids  over  a  range  of  flow  conditions  from  the 
laminar  to  the  turbulent  regime.  The  effect  of  turbulence  and  buoyancy  was  studied 
in  detail.  The  predictions  of  the  model  were  compared  with  experimental  data.  It 
was  demonstrated  that  buoyancy  may  play  a  significant  role  in  influencing  heat 
transfer  even  at  Reynolds  numbers  as  high  as  100,000.  The  direction  of  body  forces 
with  respect  to  flow  direction  was  also  found  to  be  very  important  in  determining 
heat  transfer  in  the  intermediate  range  of  Reynolds  numbers. 

Computation  of  Thermal  Stability  in  Tet  Fuels:  Advanced  thermal  stability  models 
developed  at  Wright  Laboratory  were  incorporated  into  the  code.  The  models  were 
tested  for  physical  accuracy  as  well  as  numerical  stability.  The  model  predictions 
were  compared  with  deposition  data  in  the  literature  as  well  as  data  obtained  from 
an  experimental  study  conducted  under  this  project.  The  comparisons  were  used  to 
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further  refine  and  calibrate  the  models. 


Experimental  Study:  A  detailed  experimental  study  was  performed  at  the 
University  of  Iowa  as  well  as  Wright  Laboratory  to  obtain  benchmark  validation 
data  for  the  models.  The  experiments  were  done  on  jet  fuels  (Jet-A)  as  well  as 
simulant  fluids  such  as  SF6.  The  experiments  on  fuels  yielded  data  on  temperature 
distributions  as  well  as  deposition  rates  and  the  experiments  on  SF6  yielded  data  on 

wall  temperature  distributions.  These  data  were  used  to  assess  the  validity  of  the 
model  and  to  refine  and  calibrate  model  constants. 

The  completion  of  the  Phase  II  study  has  resulted  in  a  code  that  will  enable  the 
computation  of  heat  transfer  and  fuel  thermal  stability  in  complex  aircraft  heat 
management  systems.  The  code  can  be  used  to  study  the  effects  of  various  operating 
conditions  on  the  performance  of  the  system. 

12.2  Plans  for  Future  Work 

An  immediate  extension  of  the  Phase  II  work  has  been  in  the  application  of  these 
models  to  analyze  endothermic  chemistry  in  aircraft  fuel  systems.  This  is  a  project 
that  CFDRC  is  currently  executing  for  the  US  Air  Force  in  collaboration  with 
organizations  such  as  Wright  Laboratory  and  Pratt  &  Whitney.  CFDRC  is  also 
exploring  the  potential  of  commercializing  this  code  to  organizations  developing 
aircraft  fuel  system  components  such  as  Pratt  &  Whitney,  Delavan,  AlliedSignal, 
GEAE,  etc. 

The  current  Phase  II  work  has  also  laid  the  foundation  for  developing  a  general 
purpose  capability  for  the  analysis  of  supercritical  fluids.  Many  engineering 
applications  (such  as  cleaning  solvents,  high  pressure  systems,  etc.)  involve  the  use 
of  supercritical  fluids.  CFDRC  will  pursue  these  and  other  opportunities  for 
applying  the  models  developed  during  this  study. 

12.3  Publications  Resulting  from  the  Phase  II  Study 

The  publications  resulting  from  the  Phase  II  study  are  listed  below. 
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a.  N.  Zhou  and  A.  Krishnan,  "Calibration  of  a  Global  Mechanism  for  Jet  Fuel 
Deposition,"  to  be  presented  as  AIAA-96-3238,  at  the  32nd 
AIAA/ASME/SAE/ASEE  Joint  Propulsion  Conference,  Lake  Buena  Vista,  FL, 
July  1-3, 1996 

b.  N.  Zhou  and  A.  Krishnan,  "Numerical  Simulation  of  Coupled  Heat  Transfer 
Characteristics  for  Flows  at  Supercritical  Pressure,"  to  be  presented  as  paper 
96-178,  at  the  31st  Intersociety  Energy  Conversion  Engineering  Conference, 
Washington  DC,  August  11-16,  1996. 

c.  N.  Zhou  and  A.  Krishnan,  "Laminar  and  Turbulent  Heat  Transfer  in 
Supercritical  Fluids,"  presented  at  1995  National  Heat  Transfer  Conference, 
Portland,  OR,  Aug.  5-8, 1995 

d.  N.  Zhou  and  A.  Krishnan,  "Supercritical  Heat  Transfer  in  Aircraft  Fuel 
Systems,"  presented  at  the  31st  AIAA/ASME/SAE/ASEE  Joint  Propulsion 
Conference,  San  Diego,  CA,  July  10-12,  1995. 

e.  N.  Zhou,  A.  Krishnan  and  A.J.  Przekwas,  "A  Numerical  Method  for  Reacting 
Flows  with  Multi-Step  Stiff  Chemical  Kinetics,"  presented  at  the  31st 
AIAA/ASME/SAE/ASEE  Joint  Propulsion  Conference,  San  Diego,  CA,  July 
10-12, 1995. 

f.  A.  Krishnan,  N.  Zhou  and  M.  Giridharan,  "Transport  Phenomena  in 
Supercritical  Fluids,"  presented  at  26th  AIAA  Fluid  Dynamics  Conference, 
San  Diego,  CA,  June  19-22,  1995. 

g.  N.  Zhou  and  A.  Krishnan,  "Heat  Transfer  Characteristics  of  Supercritical 
Fluids,"  presented  at  Fourth  ASME/JSME  Thermal  Engineering  Joint 
Conference,  Hawaii,  March  19-24, 1995. 

h.  T.  Edwards  and  J.  Krieger,  "The  Thermal  Stability  of  Fuels  at  480°C  (900°F). 
Effect  of  Test  Time,  Flow  Rate  and  Additives,"  presented  at  the  1995 
International  Gas  Turbine  Conference  and  Exposition,  1995. 
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12.4  Conclusions 


In  summary,  the  current  Phase  II  study  has  resulted  in  the  first  commercial  CFD 
code  that  has  advanced  models  to  simulate  flow,  heat  transfer  and  chemistry  in 
supercritical  fluid  systems.  All  of  these  models  have  been  validated  against 
benchmark  experimental  data  obtained  during  the  study. 

The  models  and  code  are  ready  to  be  applied  for  the  design  analysis  of  complex 
aircraft  fuel  systems  components.  The  application  of  these  models  will  help 
optimize  performance  of  aircraft  thermal  management  systems. 
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13. 
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APPENDIX  A 


SURFACE  DEPOSITION  DATA 


A-l 
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w 


EE 


Tube  number 


Carbon  Deposition  (pg)  per  2-inch  section 


7/94-5  7/94-6  7/94-7  7/94-8  7/94-10 


7320  2120  1570  510 


5833  1330  I  1780  614 


649  2080 


1780  832 


1930  1430 


94.3  1680 


7/94-11 


22 

97.7 

176 

104 

216 

649 

23 

106 

150 

101 

155 

385 

24 

saved 

saved 

saved 

95.9 

193 

A-4 
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Carbon  Deposition  (p,g)  per  2-inch  section 

Tube  number 

7/94-1 

7/94-1 

7/94-3 

7/94-4 

1 

13.2 

33.9 

33.8 

26 

2 

47.1 

24.4 

19.3 

20 

3 

35.6 

27.8 

16.7 

31.4 

4 

21.6 

21.4 

19.4 

30.5 

5 

24.5 

80.3 

20.6 

20.4 

6 

25.7 

309 

89.9 

17.9 

7 

23.3 

310 

150 

18.8 

8 

43 

90 

305 

17.7 

9 

155 

104 

394 

18.3 

10 

452 

217 

531 

17.2 

11 

243 

640 

126 

20.9 

12 

28.7 

140 

56.4 

17.3 

13 

31 

92.5 

42.2 

16.8 

14 

44.2 

273 

21 

18.4 

15 

62.2 

341 

26.1 

25.2 

16 

37.4 

249 

22.6 

44.5 

17 

43.8 

256 

26.9 

23.1 

18  ! 

50.1 

307 

29.7 

26.7 

19 

71.1 

26.6 

42.5 

39.7 

20 

73.5 

20.9 

68.8 

47.3 

21 

116 

23.7 

49.9 

48.8 

22 

193 

32.8 

38 

42.8 

A-5 
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APPENDIX  B 


EXIT  FILTER  DEPOSITION  DATA 


Test  Name  Filter  Deposition  (jig) 


533 


7/94-1  inlet 


outlet 


outlet 


7/94-2  inlet 


Notes 


first  5  3/4  hours 


second  4 1/4  hour 


utlet  1724 


7/94-3 


7/94-4 


There  are  no  filters  for  this  test. 


inlet 

483 

outlet 

4300 

outlet 

4460 

inlet 

374 

outlet 

4280 

outlet 

3660 

inlet 

229 

outlet 

1890 

inlet 

348 

outlet 

4820 

outlet 

2650 

inlet 

480 

outlet 

24260 

inlet 

838 

outlet 

10590 

inlet 

252 

outlet 

1280 

7/94-5 


7/94-6 


7/94-7 


7/94-8 


7/94-9 


7/94-10 


first  41/2  hours 


second  5  horns 


first  5  hours,  dirty  filter  case 


second  5  hours 


first  7  hours 


second  3  hours 


test  ended  with  pump  failure 


7  pm  filter 


7/94-11 


8/94-1 


8/94-2 


8/94-3 


8/94-4 


8/94-5 


inlet 

183 

outlet 

76.7 

inlet 

1270 

outlet 

1040 

outlet 

3030 

inlet 

538 

outlet 

460 

outlet 

saved 

inlet 

893 

outlet 

2290 

outlet 

1360 

inlet 

887 

outlet 

1790 

outlet 

1680 

inlet 

2670 

outlet 

1240 

outlet 

2730 

outlet 

2270 

first  5  hours 


second  3  1/2  hours,  power  failure 


first  2  1/2  hours 


first  6  hours 


second  4  hours 


first  7 1/2  hours 


second  5  hours,  power  failure 


first  4  hours 


second  10  hours 


next  6  hours 


APPENDIX  C 

WALL  TEMPERATURE  DATA 


C-l 
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psitjiS 


E 


m 


Wall  Temperatures  (°F) 

Distance  along  tube  (in.) 

7/94-1 

day  2 

7/94-2 

7/94-3 

7/94-4 

day  2 

5 

399 

347 

573 

384 

188 

185 

10 

814 

821 

963 

828 

538 

557 

20 

883 

918 

932 

939 

687 

30 

806 

839 

670 

673 

692 

721 

40 

818 

835 

885 

801 

673 

699 

50 

844 

858 

668 

846 

481 

472 

60 

915 

965 

895 

899 

800 

813 

70 

1002 

1140 

945 

917 

609 

815* 

80 

1055 

1119 

1121 

809 

569 

590 

92 

820 

840 

782 

921 

745 

794 

92  (fuel) 

72  7 

753 

683 

682 

485 

497 

C-2 


Distance  along  tube  (in.) 
6 


Wall  Temperatures  (°F) 


45 


979 


1083 


955 


msm 

50 

1107 

1065 

1132 

1106 

1107 

1108 

50  (fuel) 

761 

704 

738 

740 

770 

754 

■1 

Distance  along  tube  (in.) 

day  2 

8/94-3 

day  2 

8/94-4 

8/94-5 

day  3 

2.5 

397 

380 

212 

393 

355 

5 

531 

570 

535 

451 

■■ 

10 

1025 

954 

931 

940 

840 

15 

1101 

1153 

954 

1007 

20 

1375(?) 

1124 

1116 

1020 

1090 

945 

25 

972 

1212 

1210 

1015 

914 

1003 

32.5 

1071 

1033 

1035 

1115 

1042 

1150 

40 

1021 

1041 

1032 

1078 

1000 

1024 

45 

977 

1015 

1030 

940 

1015 

999 

50 

1047 

1102 

1132 

1031 

943 

946 

50  (fuel) 

765 

785 

808 

760 

800 

